Lattice Calculation of Heavy-Light Decay Constants with Two 

Flavors of Dynamical Quarks 



C. Bernard 1 , S. Datta 2 ' 3 , T. DeGrand 4 , C. DeTar 5 , Steven Gottlieb 2 , 
Urs M. Heller 6 , C. McNeile 7 , K. Orginos 8 , R. Sugar 9 , D. Toussaint 10 

(The MILC Collaboration) 



1 Washington University, 

o 

O 1 St. Louis, Missouri 63130, USA 

(N : 

o ■ 2 Indiana University, Bloomington, 

C/5 . Indiana 47405, USA 

3 Fakultat fur Physik, 
■ Universitat Bielefeld, 

: 

D 33615 Bielefeld, Germany 

O , 

\ 4 University of Colorado, 

g ; Boulder, Colorado 80309, USA 



5 University of Utah, 
Salt Lake City, Utah 84112, USA 
! 6 CSIT, The Florida State University, 



X 
S3 



Tallahassee, Florida 32306-4120, USA 
7 Dept. of Math. Sci., 
University of Liverpool, 
Liverpool, L69 3BX, UK 
8 RIKEN BNL Research Center, 
Upton, New York 11973, USA 
9 University of California, 
Santa Barbara, California 93106, USA 
10 University of Arizona, 
Tucson, Arizona 85721, USA 

(Dated: February 1, 2008) 



1 



Abstract 

We present results for f B , f Bs , /d 5 fo s and their ratios in the presence of two flavors of light 
sea quarks (Nf = 2). We use Wilson light valence quarks and Wilson and static heavy valence 
quarks; the sea quarks are simulated with staggered fermions. Additional quenched simulations 
with nonperturbatively improved clover fermions allow us to improve our control of the continuum 
extrapolation. For our central values the masses of the sea quarks are not extrapolated to the phys- 
ical u, d masses; that is, the central values are "partially quenched." A calculation using "fat-link 
clover" valence fermions is also discussed but is not included in our final results. We find, for exam- 
pie, f B = 190(7)(lf 7 )(^ 1 )(l«) MeV, f B Jf B = 1.16(1)(2)(2)(±*), f Da = 241(5)(±g)(±»)(±g) MeV, 
and /b/ Id s = 0.79(2)(1|)(3)(^q), where in each case the first error is statistical and the remaining 
three are systematic: the error within the partially quenched Nj = 2 approximation, the error due 
to the missing strange sea quark and to partial quenching, and an estimate of the effects of chiral 
logarithms at small quark mass. The last error, though quite significant in decay constant ratios, 
appears to be smaller than has been recently suggested by Kronfeld and Ryan, and Yamada. We 
emphasize, however, that as in other lattice computations to date, the lattice u, d quark masses 
are not very light and chiral log effects may not be fully under control. 

PACS numbers: PACS numbers: 12.38.Gc, 13.20.-v, 12.15.Hh 
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I. INTRODUCTION 



Accurate values for the leptonic decay constants of B and B s mesons are crucial for inter- 
preting experimental measurements of B-B mixing and bounds on, or future measurements 
of, B s -B s mixing. Knowledge of the decay constants (coupled with knowledge of the corre- 
sponding B-parameters) makes possible a determination of the CKM elements Vtd and Vt s 
from these experiments. 

In the D-meson sector, CLEO-c will measure leptonic decay rates at the 3-4% level [|T|. 
Assuming 3-generation unitarity, this translates into determinations of fu s and fr> with 
roughly 2% accuracy. Coupled with accurate theoretical computations of ratios such as 
Jb/ Jd s , this will provide crucial information about the 5-sector. In addition, if computations 
of the D and D s decay constants themselves can be performed at the few percent level, the 
experiments will directly determine V cs and V c d with similar precision. 

At least in principle, lattice QCD offers a means to compute quantities such as fs or 
Jb/Ids with control over all sources of systematic error. Here, we present a computation 
by the MILC collaboration of the decay constants fs s , /d, Id 3 , and their ratios. We 
take into account the effects of virtual quark loops from two light flavors of sea quarks; 
i.e., we have two "dynamical quarks." Additional discussion and preliminary results for 
the dynamical calculation can be found in Refs. [|2], Our earlier work, which focused on 
the quenched approximation and used dynamical configurations only for an estimate of the 
quenching errors, appeared in ||, with further details in ||. 

This paper is organized as follows: Our lattice formalism is presented in Sec. [TI|. We 
discuss the Fermilab approach to heavy quarks on the lattices || , and explain how we adapt 
it to Wilson and nonperturbatively improved clover quarks. We also explain our use of 
perturbative renormalization and the choice of scale ("g*"). 

Section |II| gives the lattice computational details. We discuss the generation of con- 
figurations, the evaluation of quark propagators for Wilson, clover, and static quarks, and 
various aspects of the analysis, including fitting and extrapolation. The most significant 
open issue here involves the effect of chiral logarithms on the light quark mass extrapola- 
tions. In important recent work, Kronfeld and Ryan || and Yamada [|9] (building on work 
of the JLQCD collaboration [K|) have argued that standard linear or quadratic extrapola- 
tions from typical lattice light quark masses miss the logarithms at low mass and drastically 
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underestimate such quantities as /b„//b- Since the quark masses available to the present 
calculation are of this typical size, it is imperative that we estimate the chiral logarithm 
effects as best we can. We devise a method to estimate, at least crudely, such effects. The 
method is based on the extrapolation of the ratio of the light-light to the heavy-light decay 
constant. 

In Sec. [TV] we revisit the quenched approximation. The dominant source of systematic 
error in our previous quenched computation f|] was the continuum extrapolation. Two new 
features of the current analysis have significantly reduced that error: (1) new running with 
both Wilson and clover quarks and (2) a new central value for the scale q* for the heavy-light 
axial current 

We then turn to the dynamical quark data in Sec. [V|. The improved control over dis- 
cretization errors in the quenched approximation gives us more confidence in the central 
value and errors deduced from the continuum extrapolation of our dynamical quark data. 
Other sources of systematic error, including the chiral extrapolation, higher order pertur- 
bation theory, and the effects of partial quenching are also discussed in detail. Finally we 
clarify the effects of "fat-link" fermions using some new test runs in the quenched case. 
These shed light on why the preliminary values for heavy-light decay constants with fat-link 
fermions on dynamical configurations were anomalously low ||12|| . 

Our conclusions and the outlook for reducing the main systematic uncertainties are dis- 
cussed in Sec. [VI]. We describe work in progress that addresses the outstanding issues in the 
chiral and continuum extrapolations. 

The computation presented in this paper is rather complicated: we use several different 
actions, operators, renormalizations, fitting techniques, and extrapolations. Part of the 
reason for this is that the simulations with dynamical quarks are extremely demanding 
computationally and therefore have taken years to complete. During that time, as our 
understanding of the physics and analysis issues grew, our methods evolved. The variety 
of methods used does have one important virtue: it allows us to estimate many of the 
systematic effects in a direct way. 
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II. FORMALISM 



In a groundbreaking paper || , El-Khadra, Kronfeld and Mackenzie (EKM) show that one 
can make sense of heavy Wilson-like 1 fermions on the lattice even when am,Q>l, where a and 
tuq are the lattice spacing and heavy quark mass, respectively. Indeed, in the nonrelativistic 
limit rriQ 3> Aqcd, they show that the effective Hamiltonian has the form (after Foldy- 
Wouthuysen-Tani transformation) 

» = 0(m 1+t ^„-^-!^)q + o(iK), (i) 

where Q is the effective heavy quark field, D is the spatial covariant derivative, B is the 
chromomagnetic field, £ are the Pauli matrices, and tuq is a generic heavy quark mass. The 
masses Mi, M2, and M3 are particular functions of the bare heavy quark mass amo that 
depend on the quark action. Here m is given by 

where kq is the heavy quark hopping parameter and k c is its critical value. 2 The "pole 
mass," Mi, controls the exponential decay of the zero-momentum propagator in Euclidean 
time, but is just an additive constant in bound state energies. The nontrivial physics of a 
heavy quark in a heavy-light bound state is controlled at this order by the "kinetic mass," 
M2, which fixes the heavy-quark energy- momentum dispersion relation, and the "magnetic 
mass," M 3 , which governs chromomagnetic effects, such as hyperfine splittings. 

For computations of heavy-light decay constants, one also needs to know how the lattice 
axial current, q'Jo'JsQ, renormalizes. At tree level but through order (1/mg), EKM show 
that the renormalization is given simply by the replacement, Q — > Qi, where the tree- level 
improved field is 



Qiix) = J2K Q e aM i [1 + adtf-D] Q(x) , (3) 



1 "Wilson-like" means that the fermion action includes the naive discretization of the Dirac equation plus 
a Wilson term to remove doublers. There may be further additional correction terms to reduce lattice 
artifacts. Standard Wilson fermions as well as "clover" fermions |l3| fall into this class. 

2 We assume throughout this paper that the spatial and temporal hopping parameters are chosen equal, 
that the Wilson parameter r = 1, and that the spatial and temporal parts of the clover term, if present, 
have equal coefficients. This is not the complete generality of Ref. M, but is sufficient for our purposes. 
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with d\ another function of am . We have included the standard ^2kq factor needed to go 
from lattice-normalized to continuum-normalized fields. 
At tree level, one has 

aMi = ln(l + amo) , 

aroo(l + amo)(2 + amo) 



aM< 



2 



aM : 



2 + AanriQ + (amo) 2 
e aMl sinh(aMi) 
1 + sinh(aMi) ' 

amo(l + am )(2 + amo) 
2(1 + am ) + cswcotio(2 + am ) 



dl 2(1 + am ) (2 + am ) ' ^ 

where c$w is the coefficient of the clover term. 

A large fraction of the one-loop corrections to eq. (|]) can be included with tadpole 
renormalization |14]]. We use k c to define the mean field value, uq, of the gauge link: uq = 



1/(8k c ). Absorbing uq into k gives a tadpole-improved hopping parameter, R, and bare 
mass, m : 

K = UqK = , K c = 1/8 

C 1 I 

arfiQ = am /u = — - — = 4— - 4 . (5) 

We denote the tadpole- improved versions of the quantities in eq. ([|) by Mi, M 2 , M 3 , and 
di. They are found simply by replacing ttiq — > ttlq. Similarly, the tadpole-improved version 
of eq. (pi) is 



Qi(x) = \/2K Q u e aM i [l + ad 1 'y-D]Q(x) 

1 _ ^9. [i + adtf.fi} Q( x ) . ( 6 ) 

For some applications, the di term here may be neglected. It is therefore convenient also to 
define 

Qf>(x) = Jl-^LQ(x). (7) 

We take the physical mass of our lattice heavy-light mesons to be the meson kinetic mass, 
Mq 9)2 . Although Mq q ^ could be determined directly from the meson dispersion relation, 
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that would require the computation of meson propagators with nonzero momenta, which in 
any case are rather noisy. Instead we define Mq^ 2 by IjUjfl 

M Qq>2 = M Qqjl + M 2 -Mi, (8) 

where Mq ?) i is the pole mass of the meson determined on the lattice, and M 2 and M\ refer 



to the heavy quark, Q. The UKQCD collaboration, in Fig. 8 of Ref. [16|, compare the 
kinetic meson mass determined by the dispersion relation with that given by eq. @ (but 
without tadpole improvement). The agreement is good, and would in fact be still better if 
the tadpole improved version were used. 

A. Wilson fermions 

For Wilson fermions (csw — 0), the magnetic mass M3 is not equal to the kinetic mass 
M 2 , even at tree level. As discussed in Refs. [0 and this produces an error at fixed a 
of 0((c mag — 1)Aqcd/ f MQ q ), where c mag = M 2 /M 3 . Hence there is little point in keeping 
the d\ term in eq. (||), which is also of 0(AQ C D/MQ q ). (Indeed, keeping such terms with- 
out including at least the 0(g 2 ) perturbative corrections to them is likely to increase the 
systematic error [[J.) We thus set d\ = in the Wilson case and use eq. (|7]). 

The 1-loop mass- dependent perturbative matching for the heavy-light axial current has 
been calculated by Kuramashi [|1^]. Since we include tadpole renormalization for both the 
heavy and light quarks through eq. (|?p (with di = 0), we adjust the result to reflect our 
choice of k c to define the mean link. (Ref. |H| uses the Feynman-gauge link.) The continuum 
contribution to the matching generates a logarithm of the heavy quark mass. In [|Tj|], this 



is taken as log Ml . Since we take M 2 as the physical mass, we also adjust the result of |18 
to replace Mi with M 2 in the logarithm. 

Additional issues for perturbative matching are the definition of the coupling constant 
and the scale at which it is evaluated. We use the coupling ay, whose value at scale 3.4018/a 



is defined in terms of the plaquette IJI4 , [l9j . It has become standard to evaluate the coupling 



at the Lepage-Mackenzie scale q* |14j], defined by 



, /A /(g)log(g 2 ^ 



l0E 'fe*» » - /A/fa) ' (9) 

where I(q) is the complete integrand for the quantity of interest. In other words, the 1-loop 
axial current renormalization constant, Za, is given by Za = 1 + av(q*)CpZA/ {^), where 



Cf is the quadratic Casimir and za = / c? 4 ? I{q)- Here, we need in three cases (light- 
light, static-light, and heavy-light), which we denote by Z^ q , Z^ tatq , and Z\^ VR , respectively, 
where kur emphasizes that we are talking about the heavy-light renormalization constant 
computed by Kuramashi We adopt corresponding notation for za and I. 

Unfortunately, at the time the analysis described here was performed, q* for Z\^ VR had 
not been determined |2fJ. For with K c -tadpole improvement, q* = 2.32/ a \2~I\, which we 
use here when fixing the scale through /„.. For Z A ts,i,q with plaquette tadpole improvement, 



Hernandez and Hill [22] found q* = 2.18/ a. Since the light-light and static-light values of 
q* were so close, it was argued in [f| that either could be used in the heavy-light case, and 
in fact the light-light value q* = 2.32/ a was chosen for the standard computation (central 
value) . 



Recently, Bernard and DeGrand [1 1 have repeated the Hernandez and Hill computation. 
They find a significantly different value of q* for the K c -tadpole- improved Z^ tatq . Their result 
depends on the heavy-light mass, which enters through the continuum part of the matching. 
However the mass dependence is rather weak over the range of masses used in the current 
numerical work, and it is therefore adequate to use an average value q* ~ 1.43/ a. Since 



Ref. |]TTJ has not yet appeared, it may be helpful to summarize here the reasons for the 
disagreement with Ref. |2"2"fl . 

First of all, Ref. sets to zero certain parts of the lattice integrand whose contributions 
to the matching vanish by contour integration. This is a standard procedure |£| for evalu- 
ating integrals involving a static quark propagator. However, such integrals do not vanish 
when the integrand is first multiplied by log(g 2 ), as in eq. @. Bernard and DeGrand argue 
that it is incorrect to discard parts of the integrand unless their contributions to both the 
numerator and denominator of eq. @ vanish. 

Secondly, there are "constant" terms in the matching coming from the dimensionally 
regularized continuum integrals. Hernandez and Hill treat these as constant over the 4- 
dimensional Brillouin zone. Similarly, the log(amg) term in Z^ tatq , which comes from both 
continuum and lattice integrals, is set to zero in p|J by the choice a = X/ttiq. In contrast, 
Ref. [DJ keeps the full continuum integrands as well as the full lattice integrands. This does 
introduce a small amount of arbitrariness: the dimensionally regulated continuum integrals 
must be replaced by finite, subtracted 4-dimensional integrals, and there is some freedom 
in how the subtraction is done. However, as long as the subtraction is "reasonable," the 
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arbitrariness in q* is small. If we accept the results of Ref. [11|, then q* for ^ tatq i s n o 
longer very close to q* for Z^. Instead, we take the static-light q* ~ 1.43/ a for Z^ UR . Since 
arrtB is quite large on our lattices (~ 1.2 to 4), we believe this is a reasonable choice. Of 
course, it is always necessary to consider a range of q* to estimate perturbative errors, and 
the range we pick (see Sec. |TV| ) includes the light-light q*, as well as the values in Ref. |20 



Summarizing the results of this section, we may express the th component of our renor- 
malized heavy-light current as 



a kuh = z^^y^QO , (10) 

with q* = 1.43/a, with Q® given by eq. (|7|), and with a corresponding expression for the 
tadpole improved light quark field, q°j. For convenience, we also use eq. ([K]) for the light- 
light pseudoscalars ("pions"), even though the mass dependence of Ref. |18| is negligible in 



that case. We take q* = 2.32/ a |2T| for light-light renormalization. 

The errors in our heavy-light Wilson calculation are formally 0(ciAqcd) and 0(ay). Note 
that, in the Fermilab formalism, one should think of these errors as multiplied by an arbitrary 
(but presumably 0(1)) function of ciMq, since we are working to all orders in aMq. Thus 
the dependence on a is in general complicated. An example of such complicated behavior 
is the difference between lattice chromomagnetic effects, which go like aAQC£>/(aM 3 ), and 
the desired behavior ciAqcd/ '(aM 2 ) (see eq. (|)). Of course, in the truly asymptotic regime 
where aMg <C 1, the leading errors are indeed linear in a, but this region is not currently 
accessible in practical calculations. 

In the static-light case, we have 

^STAT = ZY^{ q *)? llGl ,h , (11) 

where h is the static quark field, Z^ TAT is the one-loop renormalization constant for the 
static-light current p3], pi with tadpole improvement, and q* = 1.43/a. 



B. Nonperturbative Clover fermions 

For our computations with clover fermions, we take the clover coefficient csw calculated 
nonperturbatively by the ALPHA collaboration @. The th component of the renormalized, 
improved (through 0(a)) light-light axial current (which is needed here to set the scale with 



9 



f n ) is then 

A^ p = Zf y/AKg^il + b A am )[A + c A ad P 5 ] , 

A = gi7o75?2 ; P5 = qii^2 , (12) 

where K qi , n g2 are the hopping parameters of the light quarks, and m = ( m qi,o + m q 2 ,o)/2 
is their average bare mass. and c A are the nonperturbative values given in 0. The 

coefficient b A has not been determined nonperturbatively by the ALPHA collaboration, 
although the difference b A — bp has ||25|| . Bhattacharya et al. |26| have determined b A at 



ft = 6.0 and 6.2, but not at ft = 6.15, which is one of the couplings used here. Our b A is 
therefore taken from perturbation theory |p7| |, but with coupling ay(q*), with q* chosen as 



the value (= 1/a) that produces the nonperturbative result |7j for the similar quantity by. 
This gives b A = 1.47 at ft = 6.15 and 1.42 at ft — 6.0. In the systematic error analysis, we 
allow b A to vary over a range of values 3 (see Sec. |TV| ). 

For chiral extrapolations, our canonical procedure (see Sec. |I1 f"DD for both Wilson and 



clover quarks is to use the kinetic quark mass aM 2 (eq. (^)) as the independent variable. 
Through 0(a), this is equivalent to the clover "improved quark mass," am q = am (l + 
b m amo), with the choice b m = —0.5. A nonperturbative determination by the ALPHA 
collaboration |25|] gives instead b m ~ —0.709 at (3 = 6.0 and b m w —0.695 at /3 — 6.15; while 
boosted perturbation theory with the result of Ref. [^7| gives b m = —0.662 at (5 = 6.0 and 
b m —0.655 at /3 — 6.15. In the clover case, we have tried both these sets of b m values 
instead of our canonical procedure, but found only negligible changes in the central values, 
errors, and goodness of fits. For example, fs changes by less than 0.4% at f3 = 6.0 and 0.1% 
at ft — 6.15. We do not, therefore, consider the standard improved quark mass further. 

Since eq. (|T2"D is valid only through 0(a), it is likely to produce large scaling errors if ap- 
plied to the renormalization of the heavy- light axial current for heavy quarks with arriQ } o>l. 
Instead, the straightforward approach is to use the 1-loop, 0(a), perturbative matching for 
clover fermions as calculated by Ishikawa, Onogi and Yamada |28|, ^9|] . We call this approach 
"NP-IOY," where IOY refers to the authors, and NP indicates nonperturbative, because csw 
has the value given in JF], and eq. ( |T2"D is used in computing f n . 



At j3 = 6.0, reference Q gets b A = 1.28(3)(4), which is within our range. Note however, that Refs. j[> 26 
get quite different values of c A at = 6.0, indicating that effects of 0(a 2 ) and higher play a significant 
role at this coupling. 
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Since Ref. [28] uses tadpole improvement defined through k c [29|, just as we did above, 



we may summarize our application of their result as 

< P ~ IOY = ^° Y f/7o7 5 Q, + Zir^arD Q° , (13) 
where Qi and Q°j are given by eqs. (§) and (0) (and similarly for qj and g°), and where 

Z\° Y = -a v (q*)(p ( i ) +p ( i ) )/(2aM 2 ) , (14) 
with pjl\ Pa\ an( i Pa defined in p8"| . We have used the fact that the operators qi^o^^-D Qi 



and —qi^f-Djoj^Qi have equal matrix elements between zero-momentum states, as is the 
case for our evaluation of decay constants, to combine the coefficients and p^- For a 



central value of q*, we take the result from the static-light calculation of Ref. |TTJ, using the 



appropriate value of csw an d taking the heavy quark mass to be the mass of the B. This 



gives q* = 3.34/a at j3 = 6.0 and q* = 2.85/a at (3 = 6.15. Unlike pSfl , we include the d\ 



factor for the light quark in the first term of eq. (0). This is just for convenience, since d\ 
is negligible for our light quarks. We have neglected d\ for heavy and light quarks in the 
correction term (proportional to Z^ Y ) in eq. fll3|). We remark here that with our current 
data we use NP-IOY for fs, fB s but not for fo, fu s because the approximations in p8] are 
not applicable near the D mass. 

The errors in the NP-IOY calculation are formally 0(a 2 h^ CD ) and O(oy). Again, one 
should think of these errors as multiplied by an arbitrary 0(1) function of oMq, making 
the a-dependence complicated in general. For example, if Mq is held fixed and a is varied 
in the region where aMq ~ 1, linear a dependence (~ aAg CD / 'Mq) is possible. 

Although NP-IOY is a well-defined approach to heavy-light physics with nonperturba- 
tive clover fermions, it does not take advantage for heavy quarks of the nonperturbative 
information in eq. (0). An alternative, which has been used by the APE |30| and UKQCD 



|T6| collaborations, it to apply eq. (0) for moderate amg^, where it is still approximately 
justified, and then extrapolate up to the B mass. That approach has two main systematic 
errors: First, one has no guidance from HQET about the order of the polynomial in l/m,Q 
with which to extrapolate to the B mass, and second 0(a 2 ) errors, while relatively small for 
the moderate masses studied, can grow rather larger when extrapolated over a wide mass 
range. We prefer instead a different method, which we call "NP-tad" for reasons discussed 
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below. The idea here is to replace eq. ([TJ) by an equivalent expression through 0(a) but 
which has the advantage of having a reasonable limit for large arriQfi. The modified eq. flT2|) 
is then used directly at or near the B mass. 

Equation ( |T2"D applied to a heavy-light current gives 

< p = + b A amQ >° ± amq '° )[A + c A ad P 5 ] , 

A = q^oisQ ; P5 = qibQ • (15) 

Note that A^ p does not approach a static limit for amg^ — > 00. Instead, it goes to —00 
because ca < and doPs ~ sinh(MQ 9i i) ~ amQ ]0 . (^Qg,i is the meson pole mass.) Even if 
ca were zero, A^ p would still blow up because of the term ^amg^' For this reason, small 
discretization errors for moderate arriQfi may be magnified significantly if eq. ( |T5| ) is used to 
extrapolate to the B. 

To define the NP-tad alternative, we first let 

n(M Qq ) = ■ , (lb J 

{m Qfi + m qfi ){0\A \Qq) 

where Qq is a generic heavy-light pseudoscalar meson. Due to a cancellation of suiI^Mq^i) 
(from do) and the explicit mg^ in the denominator, one expects R has a finite limit as 
cf^Q,o — * 00. This is confirmed by our simulations. Then 



A^ p ' = Z2V 4 N K QV 1 + ( & a + 2c A R)am Qfi ^l + (b A + 2c A R)am qfl A (17) 

gives results for (0|Ag P '|Qg) that are identical to (0\AQ P \Qq) through 0(a). However, 
because KQarriQfl — > 1/2 as kq — > 0, eq. (|T7|) has a static limit, unlike fll2|) . 

The mass dependence of eq. HTTP is in fact very similar to the Fermilab formalism at 
tadpole- improved tree level. Indeed, from eqs. (0) and (0), the Fermilab version of eq. ([lj 
is 



^fnal = z A ^AK q K Q ^l 1 + am Qi0 / m y 1 + am,, / w A , (18) 

where Za = Z^ d uo is the renormalization constant for the axial current without the tadpole 
factor removed, and where we have dropped the d\ terms for simplicity. Hence eq. (|17D is 



equivalent to eq. (|18D but with a special (mass-dependent) value for the tadpole factor: u = 
(bA+2cAR)~ l ■ The similarity to tadpole improvement (within the context of nonperturbative 
renormalization) is the reason for the name "NP-tad." Note that, at tree level, where 
bA = Za = u = 1 and ca = 0, NP-tad is in fact identical to the Fermilab approach. 
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In practice, we put the d\ terms as well as the corresponding perturbative subtraction 
back into the axial current. Thus eq. (|T7| ) becomes 



^NP-tad = zY^/4^^1 + (b A + 2c A R)am Qfi yjl + (b A + 2c A R)am qfi A^ b (19) 
with 

Af-nb = (1 + av ( q *) p (™*>)) Ao + g lQl5 adrf.D Q-q adrf-D 7o75 Q , (20) 

where p^ ub ' ) is the IOY perturbative correction coming from the d\ terms alone, which we 
extract by comparing the complete perturbative result computed with and without the 
terms ||31|| . In eq. fl20|), d\ = di(afriQ t o) and d\ = rfi(am ?] o) in the second and third terms, 



respectively, with G?i(am ) given in eq. (f|). As before, d\{am q ^) is negligible and is included 
merely for convenience. 

Equations (|H]) and fl2"DD define the NP-tad approach. The errors of NP-tad are formally 



(9(a 2 Ag CD ) and 0(a 2 Mg). Thus, the errors could in principle be large as Mq increases 
at fixed a. The hope is that the requirement that the decay constants have a limit as 
M q — > oo has forced the a 2 Mq (and higher) terms to have small coefficients, but this is not 
proven. We emphasize that the NP-tad approach is logically neither better nor worse than 



the method of Refs. [|rj| [30| for B physics. NP-tad is an attempt to keep the nonperturbative 
0(a) information and yet include some higher effects in a smooth way, but there is no a 
priori guarantee that all or most of such effects are included. A posteriori, one can compare 
how well the methods scale with lattice spacing. The scaling of NP-tad results currently 
appears to be comparable to that seen with standard nonperturbative normalization and 



extrapolation in Ref. [16[]. 4 With only two lattice spacings in the NP-tad data and in ||T6f . 
however, this comparison is far from definitive. We therefore use the NP-IOY approach 
(whose errors are better understood), in addition to both the NP-tad method and standard 
Wilson fermions, in determining the central value and errors. 



4 An earlier version of Ref. |L6| showed considerably worse scaling, but that appears to have been associated 
more with the use of scale ro from the potential than with the normalization and extrapolation. 
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III. COMPUTATIONAL DETAILS 



A. Lattice Generation and Inversion 

Table | shows the lattice parameters used. Quenched lattices are generated using a 
standard combination of pseudo heat bath and overrelaxed [^3| updates. Successive 
configurations are separated by 200 iterations, where each iteration consists of 1 heat bath 
and 4 (9 for set H, (3 — 6.52) overrelaxed sweeps. The sets J and CP are new additions 
to the quenched lattices previously analyzed in Ref. ||. Dynamical fermion lattices were 
separated by 10 trajectories (each of unit molecular dynamics time) of the R algorithm [ pP| . 
(Set G, from HEMCGC [35j, is separated by 10 trajectories of time l/y/2 in MILC units.) 



Unimproved Wilson valence fermion propagators were generated for all these sets except 
J. (On set 5.7-large only light Wilson propagators were computed.) Because our calculations 
were initially limited by slow I/O speeds and lack of long-term storage, we performed the 
calculations of heavy-light meson propagators "on the fly" - - i.e., without storing quark 
propagators. The hopping parameter expansion of the heavy quark propagator, as proposed 
by Henty and Kenway P6| , makes this possible. Further, the expansion allows us to study 
a large number of heavy quark masses, with almost no additional expense. For this reason, 
we continued using the approach of Ref. |36| even after faster I/O and better storage became 



available. 

In the hopping-parameter approach, the light quark propagator, for a single spin-color 
source, is first computed with standard methods (red-black preconditioning; minimal resid- 
ual). The heavy quark propagator for the same spin-color source is then computed order 
by order in the heavy hopping parameter. At each order, the contribution to the meson 
propagators, summed over space, is stored to disk. The full meson propagator for any heavy 
hopping parameter, kq, can then be reconstructed after the fact by multiplying the stored 
results by appropriate powers of kq and summing over iterations as well as spin and color. 
Propagators in the static-light limit, where the heavy quark mass is taken to infinity, can 
be obtained by-product of this procedure. 

Our quark sources are Coulomb-gauge Gaussians. We run the overrelaxed gauge fixer 
until the sum of the trace of all spacelike links (normalized to 1 when all links are unit 
matrices) changes by less than 7 x 10~ 7 per pass. This takes, for example, about 435 passes 
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TABLE I: Lattice parameters. The upper group corresponds to quenched lattices; the lower group, 



to dynamical lattices with Nf = 2 staggered quarks. The set G was generated by HEMCGC [35]. 
Sets marked by * are new since Ref . [Q] . Heavy and light Wilson quark propagators were generated 
on all sets except J and 5.7-large. On the latter set, which includes lattices of various sizes, only 
light quark propagators were generated. Normal ("thin-link") clover propagators were computed 
on set J and CP1, a 199 lattice subset of CP. (See Table |n| for thin-link clover parameters.) Fatlink 
clover propagators were generated on set CPF (a 99 configuration subset of CP1) and RF (a 98 
configuration subset of set R). 



name 


(3 (am q ) 


size 


# confs. 


A 


5.7 


8 3 x 48 


200 


B 


5.7 


16 3 x 48 


100 


5.7-large 


5.7 


12 3 x 48, 16 3 x 48, 20 3 x 48, 24 3 x 48 


403, 390, 200, 184 


E 


5.85 


12 3 x 48 


100 


C 


6.0 


16 3 x 48 


100 


*CP 


6.0 


16 3 x 48 


305 


*J 


6.15 


20 3 x 64 


100 


D 


6.3 


24 3 x 80 


100 


H 


6.52 


32 3 x 100 


60 


L 


5.445 (0.025) 


16 3 x 48 


100 


N 


5.5 (0.1) 


24 3 x 64 


100 





5.5 (0.05) 


24 3 x 64 


100 


M 


5.5 (0.025) 


20 3 x 64 


199 


P 


5.5 (0.0125) 


20 3 x 64 


199 


*U 


5.6 (0.08) 


24 3 x 64 


202 


*T 


5.6 (0.04) 


24 3 x 64 


201 


*S 


5.6 (0.02) 


24 3 x 64 


202 


G 


5.6 (0.01) 


16 3 x 32 


200 


R 


5.6 (0.01) 


24 3 x 64 


200 
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on set D, and a comparable number of passes on the other sets. 

Our Wilson light quark propagators are computed for three values of K q , giving light 
quark masses (m,) in the range 0.7m s <m q <2.0m s , where m s is the strange quark mass. We 
analyze heavy-light mesons with 8 to 10 heavy quark masses per data set, with heavy-light 
pseudoscalar meson masses {Mq q ) in the range 1.25 GeV<Mgq<4.0 GeV. The heavy quark 
propagators are computed with 400 passes of the hopping parameter expansion. Figure [I] 
shows the convergence on set D of the hopping parameter expansion for heavy-light meson 
correlators at the maximum time separation (half-way across the lattice). The value of kq 
(0.1456) used in Fig. [I] gives a meson mass of Mq q ~ 1.1 GeV when the light quark K q is 
extrapolated to k u ^. Since 1.1 GeV is slightly lighter than the lightest value used in our 
analysis, we are confident that the expansion is under control. 



1.5 



^ 1.0 



% 0.5 



1 1 1 1 1 1 


1 1 1 1 1 1 1 1 1 1 1 1 1 








/ ^ — ss 




' - SL x 3000 - 


1 1 1 1 1 1 


i i i 1 i i i i 1 i i i i 



100 200 300 400 

hopping pass 



FIG. 1: Convergence of the hopping parameter expansion for heavy-light pseudoscalar mesons 
correlators on set D. kq = 0.1456; while K q = 0.1507, the lightest of the three light quarks 
analyzed on this set. The values of the smeared-smeared and smeared- local correlators at t = 40 
(half-way across the lattice) are shown with the solid and dashed lines, respectively. 

Because the heavy-light mesons must be constructed at each of the 400 orders of the 
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hopping parameter expansion, it is too expensive to sum the central point of the smeared 
sinks over the entire spatial volume, even using FFT's. Instead, we simply sum over 16 points 
in the L 3 spatial volume: the 8 points (0, 0, 0), (L/2, 0, 0), (0, L/2, 0), . . . , {L/2, L/2, 0), . . . , 
plus the 8 points obtained by adding (L/4, L/4, L/4) to each of the previous points. This 
fixes the lowest non-zero momentum which contributes to be (2, 2, 0) (and permutations) 
in units of 2 r n/L. For the heavy-light mesons studied here, these higher momentum states 
are suppressed sufficiently at asymptotic Euclidean time by their higher energy. However, 
in the largest physical volumes, sets N, O, U, and T, the higher momentum states for the 
heaviest mesons are quite close to the ground state, and we are required to go to large times 
(£min/a ~ 22-27) for the smeared-smeared propagators in order to make single exponential 
fits with good confidence levels. 

The static-light mesons have no such suppression, and the contribution of higher momen- 
tum states is limited only by their overlap with the sources. Using computed static-light 



wavefunctions |37| ) we find that the contamination in static-light decay constants from 
nonzero momentum states is small (~ 0.7%) for lattices with spatial size of ~ 1.5 fm (sets 
A, C, CP, D, E, G, H). However, on all other (larger) sets the contamination is expected 
to be large. We therefore have performed a dedicated static-light computation on those 
lattices, with relative smearing functions taken from and zero momentum intermediate 
states enforced by a complete FFT sum over spatial slices. In addition, the dedicated static 
light computation has been run on sets A and E (because the plateaus from the hopping 
method proved to be poor) and sets G and CP (as a check of the hopping method). On the 
latter sets, the two methods give consistent results, but we choose the dedicated method 
because the errors are smaller. Thus, only on sets C, D, and H is the hopping approach used 
for the final analysis of the static-light mesons. 

Standard ("thin-link") clover improved valence quark propagators have been generated 
using stabilized bi-conjugate gradient inversion |39|] (for light quarks) and the hopping pa- 
rameter expansion (for heavy quarks) on quenched sets J and "CP1," a 199- lattice subset 
of CP. In this case we have 5 light and 5 heavy quark masses, in the same range as for the 
Wilson valence quarks. For the heavy quarks we sum, on the fly, the orders in the hopping 
expansion for a given kq — i.e., this is a standard inversion, which does not allow a choice of 
arbitrary kq after the fact. However, our experience with the Wilson case leads us to believe 
that a large number of heavy quark masses is unnecessary: the behavior of fq g *J Mq q with 
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l/MQ q is quite smooth. Further, we are now able to perform an FFT sum of the meson 
propagators, so that zero momentum is enforced and contamination from excited states is 
reduced. The full Fermilab formalism allows us to choose heavy quark masses near the b 
quark mass in this case; we therefore did not need to compute static-light mesons here to 
stabilize an extrapolation. 

As explained in Sec. |11B| , we take the ALPHA collaboration J7|] values, where available, 
for the normalization and improvement constants of our clover fermions. The parameters 
used are shown in Table ||. 

TABLE II: Parameters used for nonperturbative, thin-link clover fermions (quenched configura- 
tions). 



set 


P 


#confs. 


csw 


Z A 


CA 


b A 


bA range 


CP1 


6.0 


199 


1.769 


0.7924 


-0.0828 


1.472 


1.256^1.586 


J 


6.15 


100 


1.644 


0.8050 


-0.0426 


1.423 


1.244^1.510 



B. Covariant Fits 

We need to fit correlators in time, extrapolate/interpolate in light and heavy quark 
masses, and extrapolate in lattice spacing to the continuum. In all cases except the last, the 
data is correlated, so covariant fits are preferable. As is well known however, it requires 
a large statistical sample to determine accurately the small eigenvalues of the covariance 
matrix. With limited statistics, such eigenvalues will be poorly determined and can make 
the covariant fits unstable. This is a particular problem in the current analysis because the 
large time dimension of our lattices and the fact that we fit two channels simultaneously 
means that we often make fits with 25 or more degrees of freedom. 

The technique we use to deal with the problem of small eigenvalues has many similarities 
to the methods proposed in ||(J but has some advantages in our analysis. It is based on 
a standard approach in factor analysis |4l]| . We first compute the correlation matrix (the 



covariance matrix, but normalized by the standard deviations to have l's along the diagonal) 
and find its eigenvalues and eigenvectors. We then reconstruct the correlation matrix from 
the eigenvectors, but omitting those corresponding to eigenvalues less than A cut , a cutoff. 
The resulting matrix is of course singular. It is made into an acceptable correlation matrix 
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by restoring the l's along the diagonal. Finally, the corresponding covariance matrix is 
constructed (by putting back the standard deviations), inverted, and used in the standard 
way for making the fits. 

The above technique interpolates smoothly between standard covariant fits, where no 
eigenvalues are omitted, and noncovariant (uncorrelated) fits, where all eigenvalues are omit- 
ted. Furthermore, because the correlation (as opposed to covariance) matrix is used, the 
eigenvalues are normalized, with the average eigenvalue always equal to 1. This allows us to 
make a uniform determination of which eigenvalues to keep, which is very important since 
we are dealing with thousands of fits, and it is impossible to examine each fit by hand. 
Our standard procedure is to choose A cut = 1, i.e., we drop all eigenvalues less than 1. 
The eigenvectors kept typically account for 90-95% of the total covariance. Indeed, when 
one changes how the covariance matrix is computed, (for example, by increasing the num- 
ber of configurations eliminated in the jackknife), the eigenvalues smaller than 1 generally 
change drastically with our typical sample size of ~ 100. The approach eliminates unstable, 
"pathological" fits completely. 

We have checked that the final results are not significantly affected when we keep several 
more (or several fewer) eigenvalues throughout. Furthermore, on our set with the great- 
est statistics (set CP, 305 configurations) we are able to compare with a wide variety of 
different cuts on the eigenvalue, as well as standard covariant fits where all eigenvalues 
are kept. We find that central values almost always agree within one statistical sigma, 
and usually differ by much less than that. In data discussed below (and tabulated at 
http : //www. physics .wustl . edu/~cb/Nf =2_tables) we show for comparison fits with dif- 
ferent eigenvalue cuts for this set. 



One disadvantage of the current approach, as well as to the methods in |40| , is that there 
is no true quantitative measure of "goodness of fit." When eigenvalues are removed, the 
truncated chi-squared, Xcut) tends to be (but is not always) considerably smaller than the 
"true" x 2 from a complete covariant fit. However, our experience has shown that requiring 
Xcut/d.o.f. < 1 with A cut = 1 produces fits that are almost always acceptable by a standard 
criterion (CL > 0.05, with CL the confidence level) when the data allow a fully covariant 
(uncut) fit. This is not the case for noncovariant fits (A cut = oo). Such fits may have 
Xcut/d-o.f. <C 1 and yet extend to a fully covariant fit with extremely small CL. For example, 
exponential fits to a correlator which include several points clearly outside the plateau region 
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can still have xLt/d-o.f. "C 1 when A cut = oo, but not, in our experience, when A cut = 1. 

At every stage in the analysis, we compute statistical errors by the jackknife procedure. 
The covariance matrices are also computed by jackknife. For the quenched sets, there is 
no evidence of a non-zero autocorrelation length. However, in the dynamical sets, the 
errors typically increase with the number of configurations omitted in the jackknife until 
~4 configurations are omitted. To be conservative, we determine our statistical errors and 
covariance matrices by omitting 8 configurations at a time (for both dynamical and quenched 
sets). 

C. Correlator Fits and Extraction of Decay Constants 

We compute "smeared-local" and "smeared-smeared" pseudoscalar meson propagators in 
each of three cases: heavy-light, static-light, and light-light (the last with degenerate masses 
only). These correlators are defined by 

G SL (t) = X>|<(f,t)xkO,0)|0> (21) 

X 

G S8 (t) = ^'(0|x 5 (f,t)x|(0,0)|0) , (22) 

X 

where Aq is the relevant renormalized current, namely A^ VR (eq. (|10|) ), v4q TAT (eq. (|TTD ), 
a np-ioy ^ ^NP-tad ^ or ^np ( eq _ for light _ light quan tities only). 

For Gss, the prime on the sum indicates that only N sin ^ points on a time-slice are included. 
As discussed above, for all the Wilson heavy-light and light-light data, iV sink = 16. For the 
clover and the dedicated static computations, the complete sum is performed (with FFT), 
so one has A^ sin k — V = n x n y n z . In eqs. ( pT|) and (^2|), X5 is the Gaussian pseudoscalar 
source, given by 

X5 (x,t) = J2e~ iP/r "e^^q(x + y) l5 Q(x + z) . (23) 

yz 

The width r varies from 2.33 lattice spacings (set A) to 8 spacings (set H) and is chosen 
to be roughly 0.35 fm. For the Henty-Kenway hopping calculations, the sums in eq. ( |23|) 
run over even y, z only, so that we may exploit an even-odd decomposition. In eq. ( |23|) 
and below, we use the notation of the heavy-light case (quarks Q and q) generically: for 
light-light formulae, let Q —> q\ for static-light formulae, let Q — > h. 
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For large Euclidean time t, Gsl and Gss are fit simultaneously and covariantly to single 
exponential forms, with the same mass in both channels 



G 



SL 



CsLe 



-Mt 



G 



ss 



Csse 



-Mt 



(24) 



In other words, these are fits with three-parameters: M, (sl and (ss- Central values use 
A cut = 1.0 throughout, except for set CP, where A cut = 0.1. Typical effective mass plots 
for the light-light and heavy-light cases are shown in Figs. |2], |3| and [|, respectively. Here 
and below, we generally choose Nf = 2 and clover data for the plots because the quenched 
Wilson data has been discussed in more detail previously 0, |5j . 
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FIG. 2: Light-light effective masses for set R, k = 0.159. The fit ranges are 8 to 31 for Gsl and 
10 to 31 for Gss- The smeared-smeared masses are shifted upward for clarity. The long horizontal 
lines show the fit value of the mass. The error in the fit value is indicated at the left end of the fit 
lines by two short horizontal lines. 
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FIG. 3: Same as Fig. |2|, but for heavy-light effective masses on set J, kq = 0.098, K q = 0.1347. 
The fit ranges are 12 to 18 for Gsl and 8 to 16 for Gss- 

We vary the fit range (in t) in each channel over several choices that appear to be in 
the asymptotic, "plateau," region for the effective mass. Combining such choices for the 
light-light, heavy-light and static-light cases, we have approximately 25 different versions of 
the analysis on each data set. Our central values are taken from the version which has the 
best blend of small statistical errors and low, or at least acceptable, values of Xc U t/d-°-f-- 
Here "acceptable" is defined, with few exceptions, as Xcut/d-o.f. < 1-3, with A cut = 1. (In 
about 85% of the central value choices, Xcut/d-o.f. < !•) The exceptions are: (1) a few fits 
to heavy-light mesons outside the ranges of mass that we include in our final determinations 
of the decay constants, (2) a few fits to heavy-light mesons on the largest lattices, where the 
slow approach to the asymptotic regime for smeared-smeared correlators computed with the 
Henty-Kenway approach (see Sec. |111 A| ) left us with somewhat noisy data at large t, and (3) 
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FIG. 4: Same as Fig. [2], but for heavy-light effective masses for set M, kq = 0.120, K q = 0.160. 
Fit ranges are 13 to 31 for Gsl and 19 to 31 for G,ss- The heavy quark is computed with the 
Henty-Kenway hopping expansion. The late plateau in Gss is due to contamination from non-zero 
momentum states because the sink is not summed over complete time slices. 

the fit to the static-light meson with heaviest light quark on set G, which had rather noisy 
data. In these exceptional cases, we relax our definition of acceptable Xc Ut /d.o.f. to be less 
than 1.9, 1.7, and 1.6, respectively. 

Of course, with the approximately 900 channels we fit, one should expect that some fits 
over truly asymptotic ranges will have poor x 2 /d.o.f. simply because of statistical fluctua- 
tions. However, since our Xcut generally underestimates the full x 2 , we have tried to make 
choices which are more conservative than a standard criterion, of say, confidence level > 0.05. 

From the fits, the pseudoscalar decay constants, /q 9 for given quark masses are then 
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obtained via 



f Q ,VW, = 4J^^=- (25) 
as V v yi^ss 

where we use the definition of the decay constant that gives f n = 130.7 MeV: 

<0 |A c ont | Qq,P = 0) = -if Qq M Qq , (26) 

with A cont the continuum axial current. 5 

Data for the masses and decay constants for each of the sets listed in Table | are posted 
at http : //www. physics .wustl . edu/~cb/Nf =2_tables . The files "name_qq.dat," where 
"name" is the set name (A, B, E, . . . ) give dimensionless light-light masses and decay con- 
stants as a function of hopping parameter. Similarly "name_Qq.dat" and "name_Statq.dat" 
give dimensionless masses and values for a 3 / 2 /q 9 ^Mq q for heavy-light and static-light 
mesons, respectively. For set CP, additional files, with "lambda-cut=X" appended to the 
name, show the effect of various truncations of the correlation matrix. For heavy-light 
mesons the masses tabulated are the shifted masses aMQ q> 2, eq. (||); for the static-light case 
they are simply the pole masses. Time ranges, number of degrees of freedom, and Xcut are 
included for all fits. For sets CP1 and J, the nonperturbative clover lattices, "name" gets a 
further qualifier in the heavy-light case, which is either "NP-IOY" or "NP-tad" for the two 
types of renormalization performed. 

To enable the reader to see the effects of various renormalizations used, as well as to make 
possible reanalysis of the data by other groups, we have tabulated additional raw data. For 
all sets, we have separately computed correlators of the bare lattice axial current A Q : 

G b s T(t) = 5>|4>(2,f)x!(ff,0)|0) ; A = q lol5 Q. (27) 



We fit G b g[ e simultaneously with Gss, as in eq. (|24[) , giving us the three quantities M, ($f e 
and (ss- We then define S bare by 



/ IV /^bare 



v VC, 



ss 



5 When eq. ( p5f ) is used for static-light mesons with the Henty-Kenway hopping approach, an extra factor 
of ^/2 is required on the right hand side. This arises from the fact that the highest momentum state 
(tt, 7t, 7r) aliases the zero momentum state with our even-site-only source, and the higher momenta are not 
suppressed in Euclidean time for the static case. 
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S bare is basically the unrenormalized decay constant. For example, in the case of our heavy- 
light Wilson data, for which the renormalized current is given by eq. (|10]), we have 

a 2 V V 

The masses and quantities S bare are posted in the files "netme_qq_bare.dat," 
"name_Qq_bare.dat," "name_Statq_bare.dat" . 6 

For heavy-light mesons with the nonperturbative clover action, we have posted additional 
intermediate data. We define 

< m4 = q l0l5 a^-DQ (30) 
Aft ee A + (di(am Q>0 ) + ^(am.o))^ 4 (31) 
A d ^ b = (1 + a v (q*)p [ r h) )Ao + (di(omg l0 ) + dx(am ? , )) A* m4 (32) 



A™ p ee Jl + (b A + 2c A R)am Qfi Jl + (b A + 2c A R)am qfi A d ^ ub , (33) 



where "imp" stands for "improved," and where eq. (|32[) is a rewriting of eq. (|20[) , using the 
fact that the operators qi'jojs'y-D Qj and — qi'y-DjojsQi have equal zero-momentum matrix 
elements. For each current in eqs. ([301) -([33]), we define a corresponding S, as in eq. (p8|) . Re- 
sults for S dim4 , ~ dl ^ ub and S imp are tabulated in the files "name_Qq_intermediate.dat", 
where "name" now is just CP1 or J. From eqs. (|13|), (§), ©) an d (0), these quantities are 
related to the decay constants in the NP-IOY case by 



y/2 [ 3k q f 3^ 



7 IOY | 7 IOY^dim4 

Z/ 4 ^ -f- Z/ 12 — 



(34) 



Similarly, from eqs. (|19|), fl32|) and ([33"D, we have in the NP-tad case: 

V2 



^M^ q = \ ZY y^T q S im P (35) 
V2 



CL2 



3 

a 2 



Zfv^^Jl + (6 A + 2c A R)am Qfl Jl + (6 A + 2c A R)am qfi E d ^ ub . (36) 



In practice, when eqs. (|4j) and fl25"D are used to compute fQ q ^Mg q , eqs. ([34]) and (|36|) are 
obeyed only up to small corrections. This is because two separate fits are performed to com- 
pute the two terms in eq. (|34]) ; whereas the quantities are added first and then fit in eq. 



6 Again, for static-light correlators computed with the hopping approach, we include an extra factor of y/2 
on the right-hand side of eq. (|28|). In this case the heavy- light masses tabulated are the (pole) masses 
directly from the fits, not the shifted masses. 
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For eq. (|36|) , the discrepancy is due to the fact that the factors like + (pA + 2cAR)aiiiQ :0 
and their errors are not included in the fit here, but are factored in later. 

Finally, we also compute smeared-local light-light vector meson propagators. These are 
fit covariantly to single exponentials (two parameter fits). Raw data for the vector channel 
appears in the files "name_qq-vector.dat". 



D. Chiral Extrapolations 

Chiral extrapolations/interpolations are needed for the light-light pseudoscalars, which 
are used to set the scale (through f n ) and to find k c and the physical values of n U) d and k s , 
the hopping parameters of the up/down and strange quarks. (We generally determine k s 



by adjusting the degenerate pseudoscalar mass to ^2m 2 K — m*., the tree- level chiral pertur- 
bation theory value.) The light-light vectors provide alternative determinations of the scale 
(through m p ) and k s (through m^) and require additional chiral extrapolations. The heavy- 
light and static-light masses and decay constants also need to be extrapolated/interpolated 
in light quark mass to the up/down and strange quark masses. 

We have tried chiral extrapolations using either the bare light quark mass m q0 or the light 
quark tadpole- improved kinetic mass M q ^ as the independent variable. 7 For both Wilson 
and clover quarks, the confidence level of linear fits to Mq q , fg q , and f qq are better with 
M qt 2 than m qQ , so we use it from now on in all cases. 

An important question is which functional form one should fit to. Unfortunately, as in 
other lattice computations to date, we have been forced to work at fairly large values of light 
quark mass. In this region, our data for decay constants, both fg q , and f qq , are quite linear. 
There is little evidence for chiral logarithms, which should introduce significant curvature 
as one approaches the chiral limit, as emphasized recently by Kronfeld and Ryan || and 
Yamada M. This is presumably not a problem with chiral perturbation theory (xPT), 
but simply an indication that higher order terms (e.g., terms quadratic in quark mass) are 
as important as the chiral log terms in the current mass regime. Further, chiral log fits 
would introduce yet another parameter in the heavy-light case, the B-B*-tt coupling g. 8 It 



7 As discussed in Sec. II B, we also tried the standard 0(a) improved quark mass in the clover case, but 
the fits were not significantly different from those with M Qj 2- 

8 CLEO |42] has recently measured the D* + width, which gives, using lowest order %PT, a D-D*-tt coupling 
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therefore seems clear to us that fits of /q 9 to the NLO %PT form would require at least 
four parameters: the value in the chiral limit, a linear slope, the coefficient g 2 , and a higher 
order (quadratic?) term introduced to cancel most of the curvature of the logarithms in our 
relatively high mass region. With only three light quark masses on most sets, it is clear 
that such an approach is not feasible at present. In work in progress [4~5|, however, we are 
consistently using five light quark masses and expect that we will be able to include chiral 
logs and quadratic terms directly in the fits for central values. 

For our "standard" chiral extrapolations, we thus consider only quadratic and linear fits 
in M Qi 2- For each physical quantity, we choose one of these fits for the central value, and the 
other is taken, where appropriate, to give a standard chiral extrapolation error. Note that 
have similar chiral log effects in full QCD: 

1 



f B |H and U 



1 



16n 2 f 
1 

16vr 2 / 2 



, 3(l + 3ff 2 ) 2 2 , 
> n ( 4 m , M™ J + • • • ) 



(-2m 2 ln(m 2 ) + ---) 



(37) 
(38) 



Since we fix the lattice scale with f n and always use the same type of chiral fit for both fg q , 
and f qq , it is not unreasonable to expect that much of the systematic effect coming from not 
including the curvature of eqs. ( |3~7"D and fl38|) will cancel. In Sees. |IIIE| and |VD| , we explain 
how we test this assumption and estimate the chiral logarithm effects in the dynamical case. 

Figure [5] shows the chiral extrapolation of m 2 to find k c , with both linear and quadratic 
fits. Since the independent variable, M qt2 , itself depends on k c , such fits have to be iterated 
two or three times to find a self-consistent value of k c where both m 2 and M q ^ vanish. This 
has been done only for the quadratic fit in Fig. |5] to emphasize the difference with the linear 
fit. 

Table [TTT] shows the results for k c and xlut values for all the data sets. Note that the 
linear fits are uniformly very poor; while quadratic fits are quite good wherever there are 
enough light masses to compute a xLt- Although using am as the independent variable 
actually reduces the xlut va hies of the linear fits somewhat, they remain very poor. These 
features agree with what was found in Ref. H, where we performed the light-light analysis 



0.35. However, NLO %PT and the D* -> Dtt decay gives g 2 « 0.07 |43 ; while NLO %PT on the 



D* + width (eq. (21) in Q) gives g 2 w 0.22. A recent lattice computation gives g 2 « 0.45. There is 
also some uncertainty in going from the D to the B system. 
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FIG. 5: (am,) 2 vs. aM q ^ 2 for set CP1. The solid line is a quadratic fit; the dotted line, the 
alternative linear fit. The fits have x^ ut /d.o.i. = 0.06 and 36.5, respectively, with A cu t = 1. The 
bursts are the extrapolated points where (am^) 2 = 0. 

on a large quenched data set at (3 — 5.7, using five light quark masses. We therefore use 
only quadratic fits/solves from here on for m 2 vs. aM q p- Table |IV| gives the values of n s 
resulting from these fits. 

The case of f qq vs. aM q ^ is a more difficult one. Figures || and |7|, and Table [V| show the 
extrapolations. Although the linear fits in both figures appear reasonable to the eye, that 
in Fig. [7| has a rather high value of X™t/ d.o.f , as do many of the other linear fits in Table [V]. 
Where comparisons can be made, the quadratic fits are better. This was also true of the 
large quenched data set at (3 = 5.7 analyzed in Ref. [|J. On the other hand, the quadratic 
fits often have quite large statistical errors, especially in the Wilson case where there are 
only three light quark masses. Furthermore, on the sets with the finest lattice spacings (C, 
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TABLE III: Extrapolations of a 2 m^ vs. aM q ^ to find k c . When the number of degrees of freedom 
is (quadratic fit with 3 k values), a solver is used instead of a fitter. For fits, the cutoff on 
correlation matrix eigenvalues (A cu t) is 1. For set 5.7-large, which consists of several lattice sizes, 
the central value is the average over all sets, the error is a combined weighted error, and the x 2 an d 
d.o.f. shown are the ones from the volume 20 3 x 48. See Sec. VB for a description of the fat-link 
computations. 



name 




Xcut 


d.o.f. 




Xcut 


d.o.f. 


quadratic fit 


linear fit 


quenched Wilson 


A 


0.169433(237) 







0.168607(89) 


29.3 


1 


B 


0.169340(100) 







0.168383(52) 


227.0 


1 


5.7-large 


0.169748(24) 


0.3 


2 


0.168862(33) 


703.3 


3 


E 


0.161397(124) 







0.161046(89) 


20.1 


1 


C 


0.157228(95) 







0.156778(46) 


62.5 


1 


CP 


0.157274(74) 







0.156906(25) 


44.0 


1 


D 


0.151825(55) 







0.151663(35) 


28.0 


1 


H 


0.149368(20) 







0.149248(15) 


18.9 


1 




N f = 


2 Wilson 






L 


0.169422(61) 







0.168160(50) 


299.4 


1 


N 


0.169515(56) 







0.168513(25) 


496.4 


1 





0.167197(41) 







0.166483(30) 


306.8 


1 


M 


0.165919(59) 







0.165211(33) 


192.3 


1 


P 


0.165257(47) 







0.164864(28) 


144.5 


1 


U 


0.163065(27) 







0.162570(16) 


426.6 


1 


T 


0.161528(22) 







0.161887(27) 


155.9 


1 


S 


0.161400(20) 







0.160802(12) 


447.2 


1 


G 


0.161158(72) 







0.160821(45) 


38.7 


1 


R 


0.161167(23) 







0.160798(12) 


328.9 


1 


quenched clover 


CP1 


0.135342(18) 


0.1 


2 


0.135168(8) 


109.4 


3 


J 


0.135862(20) 


0.0 


2 


0.135792(16) 


45.8 


3 


fat-link clover (N = 10, c = 0.045) 


CPF 


0.125558(22) 


0.2 


2 








RF 


0.125666(25) 
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TABLE IV: Values of k s , the hopping parameter of the strange quark, from fits to the light-light 



m 2 ; in 



pseudoscalars and vectors. In the former case, we adjust the pseudoscalar mass to ^2m? K 
the latter, we adjust the vector mass to m^. The values of x 2 and degrees of freedom can be found 
by referring to the corresponding fits in Tables 111 and Data from set 5.7-large is combined as 
in Table ID. See Sec. [V B| for a description of the fat-link computations. 



name 


k s (pseudoscalars) 


k s (vectors) 




quadratic fit 


linear fit 


quenched Wilson 


A 


0.164331(432) 


0.163173(407) 


B 


0.163629(355) 


0.163709(215) 


5.7-large 


0.163916(100) 


0.163456(78) 


E 


0.158203(170) 


0.157729(351) 


C 


0.154567(97) 


0.154780(229) 


CP 


0.154857(99) 


0.154638(152) 


D 


0.150395(66) 


0.150316(102) 


H 


0.148415(42) 


0.148322(87) 


N f = 2 Wilson 


L 


0.164114(201) 


0.164064(184) 


N 


0.164436(213) 


0.164027(202) 





0.162938(121) 


0.162691(120) 


M 


0.161778(205) 


0.161830(206) 


P 


0.161518(141) 


0.161374(106) 


U 


0.159690(81) 


0.159373(103) 


T 


0.158610(121) 


0.157243(148) 


S 


0.158633(68) 


0.158423(50) 


G 


0.158795(94) 


0.158350(146) 


R 


0.158736(55) 


0.158465(85) 


quenched clover 


CP1 


0.133882(41) 


0.133515(63) 


J 


0.134665(45) 


0.134400(56) 




fat-link clover (N = 10, 


c = 0.045) 


CPF 


0.123206(46) 


0.123481(137) 


RF 


0.123440(62) 


0.123236(162) 
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CP, D, H, G, R, CP1, J) the linear fits are generally acceptable (set R is an exception). 
For these reasons, we use the linear fits for the central values and take the quadratic fits 
to estimate the "standard chiral systematic error." As mentioned above, the fits in future 
work [J5| will include chiral logs as well as quadratic terms. 

0.09 i — i — i — i — i — I — i — i — i — i — I — i — i — i — i — I — i — i — i — i — I — i — i — i — r 




04 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 L 

-0.02 0.00 0.02 0.04 0.06 0.08 



aM 



q>2 



FIG. 6: The light-light pseudo-scalar decay constant, f qq , vs. aM qt 2 for set J. The solid line is a 
linear fit; the dotted line, the alternative quadratic fit. The fits have Xcut/d-o.f. = 1.1 and 0.005, 
respectively, with A cu t = 1 (1 eigenvector of 5 kept). The bursts are the extrapolated points where 
Mg t 2 takes its physical value (i.e., n = n u ,d)- The burst on the dotted line has been displaced 
slightly for clarity. 

Figure || shows a typical linear extrapolation of the lightdight vector mass to the physical 
point for up/down quarks. As seen in Table 0, the linear fits are almost always quite good, 
and we use the scale set with m p in this way as an alternative to that from /„-. The same 



fits also give an alternative value for k s , using m<j>. The results are shown in Table [TV 
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TABLE V: Extrapolations of af qq and m p vs. aM q i to find the scale, a" 1 , using f n = 0.1307 GeV 
and m p = 0.768 GeV. When the number of degrees of freedom is (quadratic fit with 3 k values), 
a solver is used instead of a fitter. For fits, the cutoff on correlation matrix eigenvalues (A cu t) is 1. 



Data from set 5.7-large is combined as in Table III. See Sec. VB for a description of the fat-link 
computations. 



name 


a' 1 (GeV) 


Xcut 


d.O.f. 


a' 1 (GeV) 


Xcut 


d.O.f. 


a- 1 (GeV) 


Xcut 


d.O.f. 


f qq , linear fit 


f qq , quadratic fit 


m p , linear fit 



quenched Wilson 



A 


1.391(66) 


5.3 


1 


1.586(122) 







1.413(32) 


1.5 


1 


B 


1.311(48) 


6.1 


1 


1.453(66) 







1.488(25) 


0.0 


1 


5.7-large 


1.339(14) 


3.3 


3 


1.388(17) 


0.3 


2 


1.450(8) 


9.9 


3 


E 


1.780(43) 


3.4 


1 


1.851(63) 







1.848(73) 


0.1 


1 


C 


2.124(54) 


1.8 


1 


2.266(96) 







2.414(102) 


0.5 


1 


CP 


2.211(48) 


1.2 


1 


2.339(124) 







2.333(51) 


0.2 


1 


D 


3.151(91) 


2.8 


1 


3.331(111) 







3.389(108) 


0.2 


1 


H 


4.388(121) 


0.5 


1 


4.490(146) 







4.489(150) 


0.0 


1 


N f = 2 Wilson 


L 


1.375(30) 


15.7 


1 


1.585(43) 







1.545(23) 


0.0 


1 


N 


1.432(35) 


2.0 


1 


1.524(53) 







1.524(25) 


0.3 


1 





1.568(25) 


2.3 


1 


1.654(55) 







1.685(20) 


0.9 


1 


M 


1.608(47) 


1.0 


1 


1.680(56) 







1.789(37) 


0.3 


1 


P 


1.713(38) 


2.0 


1 


1.773(50) 







1.853(22) 


2.3 


1 


U 


1.800(25) 


3.2 


1 


1.868(34) 







1.888(25) 


0.3 


1 


T 


1.800(37) 


24.3 


1 


1.939(40) 







1.839(20) 


0.7 


1 


S 


2.038(28) 


7.6 


1 


2.192(31) 







2.157(19) 


0.5 


1 


G 


2.243(38) 


4.3 


1 


2.377(63) 







2.242(52) 


0.0 


1 


R 


2.194(26) 


8.7 


1 


2.306(39) 







2.269(31) 


3.6 


1 


quenched clover 


CP1 


1.994(29) 


8.0 


3 


2.093(65) 


3.9 


2 


1.908(26) 


0.2 


3 


J 


2.447(45) 


3.3 


3 


2.545(90) 


0.0 


2 


2.383(39) 


1.4 


3 






fat-link clover (iV = 10, c = 


0.045 


) 


CPF 


1.847(18) 


14.4 


3 


1.923(28) 


0.5 


2 


2.114(61) 


6.1 


3 


RF 


1.873(27) 


1.4 


1 


1.939(69) 







1.920(47) 


0.3 


1 
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FIG. 7: Same as Fig. |6|, but for set R. The linear fit (solid line) has Xcut/d-o.f. = 8.7, A cut = 1 (1 
eigenvector of 3 kept); while the quadratic fit (dotted line) has no degrees of freedom. 

Sample chiral fits for heavy-light masses and decay constants are shown in Figs |9| 
and [It]. Although the extrapolated data for all the sets are too extensive to tabulate 
here, they are available at http : //www. physics .wustl . edu/~cb/Nf=2_tables, in the files 
"name_chiral_mass.dat" and "name_chiral_fsqrtM.dat." In the latter files we give /q 9 a/ Mq q , 
rather than fQ q , because /q 9 a/ Mq q is what we will need later to extrapolate/interpolate to 
the mass of the B and D mesons. For the decay constants from sets CP1 and J, "name" 
includes the qualification NP-IOY or NP-tad because /q q will of course depend on how the 
renormalization is done. 

For almost all data sets, the linear chiral fits of the heavy-light masses are quite good, at 
least in the important range of meson masses between the D and the B. Sets L and T are 
exceptions, which is perhaps not surprising since their correlators are quite noisy to begin 
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FIG. 8: Linear extrapolation of the light-light vector meson mass to the physical point (p meson, 
indicated by the burst). The data is from set S. The fit has x^ ut /d.o.f. = 0.5 with A cut = 1. 

with, making it difficult to find good plateaus. Indeed, for set T the data is noisy enough 
that covariant chiral fits for heavy-light masses did not converge, with any choice for A cut 
except A cut = oo, i.e., noncovariant fits. However, since the linear fits to heavy-light masses 
were fine on the vast majority of the sets, we believe it is reasonable to use them exclusively. 
Linear fits to the static-light masses are always acceptable. 

The situation for the heavy-light decay constants is similar to that for the light-light decay 
constants. Again, there is a small amount of curvature in these plots, and the direction of 
curvature is the same as for the light-light case. (Compare Fig. [l(| with Figs. |] and 
Therefore, quadratic chiral fits of /q 9 are better than linear ones where the comparison 
can be made, but, as before, the quadratic fits/solves lead to significantly larger statistical 
errors. The main difference with the light-light case is that the linear fits typically improve 
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FIG. 9: Linear chiral extrapolation of the heavy-light pseudoscalar meson mass, to the physical 
point (k = K Ut d, indicated by the burst). The data is from set U, kq = 0.113, The fit has 
Xcut/d-o.f. = 0.7 with A cut = 1. 

as the mass of the heavy quark increases, so that by the time the physical b quark mass is 
reached they are generally quite reasonable. The only exceptions to this rule are the sets 
CP, for which linear fits are poor over the whole heavy-quark mass range, and T, which 
is noisy and again requires noncovariant fits. In the static-light case, linear fits are always 
good. We thus choose linear fits everywhere for the central values but use quadratic fits 
for the heavy-lights in estimating the systematic error in the standard chiral extrapolation. 
As discussed above, we expect that extrapolating both light-light and heavy-light fits in 
the same way will cancel at least some of the systematic error associate with curvature and 
chiral logarithms. Therefore, we change from linear to quadratic fits in both cases at once 
when we make our estimate of the "standard" chiral error. 
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FIG. 10: Linear (solid) and quadratic (dotted) chiral extrapolations of /g g to the physical point 
(k = indicated by the bursts). The data is from set CP1, kq = 0.09, with NP-IOY renormal- 
ization. The linear fit has Xcut/d-o.f. = 1-0; the quadratic, 0.3. The parameter A cu t = 1 in both 
cases. The burst on the dotted line has been displaced slightly for clarity. 

To summarize: For central values, we use quadratic chiral fits in aM q ^ for m\ and 
linear fits for heavy-light and static-light masses and all the decay constants. We call this 
combination "chiral choice I." Our standard chiral systematic errors are found by comparing 
the central values with the results of "chiral choice II:" quadratic fits for m\ and both light- 
light and heavy-light decay constants, and linear fits for heavy-light masses and static-light 
masses and decay constants. For the light-light vector meson masses, which enter only in 
various systematic error estimates (alternative scale determination from m p , alternative k s 
determination from m^), we always use linear fits. 
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E. Chiral Logarithm Effects 



The standard chiral systematic error just described does not directly take into account 
the sharp curvature in decay constants at very small quark mass caused by terms of the form 
— m^ln(m^) in eqs. fl3~7]) and (|38|). Putting aside the issue of scale choice, an extrapolation 
in the full theory that ignores the chiral log in fQ q is expected to overestimate fs and 
underestimate fs s //b | (since — m 2 In m 2 is concave down with rapid variation at small 
mass). However, because we set the scale with and use the same extrapolations for f qq 
and fQ q , the effect on the individual decay constants is less clear. Even ratios such as fsjfs 
are affected indirectly by the scale choice, through the fixing of k s1 the hopping parameter 
for the strange quark. (Fixing Hf,, the bottom quark hopping parameter, has little effect 
on the ratio but does represent another scale effect on the individual decay constants.) It 
is easy to see that our scale choice should push Jb and /b 3 //b in the opposite direction 
of the fg q extrapolation. We can thus hope that such effects largely cancel, and this is a 
justification for taking our central values and errors from the standard linear and quadratic 
chiral fits described in the previous section. 

However, to test the above assumption and estimate the errors induced by not directly 
fitting with chiral log forms, we need alternative methods of evaluation that do not involve 
chiral extrapolations of individual decay constants. One approach that takes advantage of 
the fact that the chiral logs in and f n are similar in magnitude is to perform chiral 
extrapolations only on the ratio /#/ f n (more precisely, fq q / f qq ). Given the expected range 
for the parameter g 2 , the ratio has a chiral log term of opposite sign from that of /b alone 



and either comparable or greatly reduced magnitude (see eqs. (|37|) and (|38|) and the footnote 
shortly before). In practice, since the slope of f gq is greater than fQ q , f qq / fQ q is more linear 
than fQ q / f qq , and we work with the former. 

Fig. |TT] shows a chiral extrapolation of f qq /fQ q for set R, with kq chosen so that Mg q is 
near the B mass. There is clear curvature, so a linear fit is not appropriate, and we fit (solve) 
quadratically. To the extent that a residual chiral log remains in the ratio, the quadratic fit 
should somewhat overestimate f n / fs and hence also /b s //b- 

However, to take real advantage of the presumed reduction of chiral logs in /„./ fs and 
known sign of the error in /b s //b, we must eliminate the dependence of the scale and k s 
on a chiral extrapolation. This means that common scale choices such as /„., m p or 
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FIG. 11: Quadratic chiral extrapolation of f qq / fg q from set R, with kq = 0.125. The bursts show 



the extrapolation to k s and K Ut d, with these values determined by method 6 in Table [VI 



cannot be used. Further, we are reluctant to employ static potential quantities such as the 
string tension, r or r x because their physical values are only known phenomenologically, 
with uncertain errors. Instead, we look at three more-or-less physical quantities associated 
with the s quark: the vector meson mass m^, and the mass and decay constant of a would-be 
ss pseudoscalar, m ss and f ss . To be precise, the ss meson is made of two valence quarks 
with the physical strange quark mass, but in a standard sea quark background: either the 
physical Nf = 3 case (with m u = rrid = m phys and mf ys ), or, corresponding more closely 
to our simulations, the Nf = 2 case (with m u = = m phys ). Here the superscript "phys" 
stands for the physical mass, and we neglect isospin violations as usual. 

The quantity / ss can be related to fx and f n in 1-loop (NLO) partially quenched chi- 
ral perturbation theory (PQ^PT), in a manner independent of the analytic p 4 ("Gasser- 
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Leutwyler") constants. Using the formulas in Refs. j48| and [g9| for Nf = 2 and Nf = 3, 
respectively, we find: 

x( 2 ) f i 

= 1 + 16Vf ln ( m -/ m ' ^ ~ 5 m « + X) « 0-93 (39) 

f (3) /• i 

^ = 1 + (t™? InK/A 2 ) - InK/A 2 ) - ln(mL/A 2 )) « 0.95 (40) 

where, in the numerical evaluation, we have used / = f w = 130.7 MeV and m 2 ss = \m 2 r] — \ml. 
(The latter relation is theoretically convenient here because it makes eq. ( flO|) explicitly 
independent of the chiral scale A.) 

For m ss , we consider a few choices. Two tree level relations have already been mentioned: 



m 



I = 2m 2 K - m\ (41) 



™ 2 SS = \ml - \m\ (42) 



With Refs. |48|, |49fl , we can also derive 1-loop relations for m ss similar to eqs. ( p9|) and 
(flOl), although they do involve analytic terms indirectly through the quark mass ratios. In 
particular, we use the Nf = 2 result: 

= ^ m /^ a [l + - J— (-ml Hm 2 Jml) + m 2 ss - ml) 1 « 0.17 , (43) 

where we have taken m s /m = 24.4 ||50|| . 

We then perform a series of analyses. For each, we choose 2 of the quantities f ss , m ss , 
and m^, and a method of evaluation for the "physical" quantities f ss (either eq. ( j39|) or 
eq. (|4"0D) and/or m ss (either eq. (f4T|), eq. ( f42|) or eq. (f£3|)). We then fit the ratio of the 
2 chosen quantities as a function of light quark mass (aM q ^) in order to determine k s . 



Generally, only an interpolation or mild extrapolation is required. Figures [L^ and [13] show 
two examples of such fits, for m ss /f ss and m ss /m^, respectively. The former uses eqs. ( f4"T|) 
and fl3"5p; while the latter, eq. (fJ^). Note that a slight extrapolation is required to find k s 
in Fig. |12[ In Fig. [13|, k s is determined by an interpolation, which is more similar to the 
standard analysis. The difference between the two situations gives some indication of the 
errors of the procedure. In Tables [VI| we show the results for k s for 12 different versions 
of such fits, performed on two different dynamical sets of configurations (the ones with the 
lightest sea quark masses). 

Given k s , the next step is to determine the scale, a -1 . We consider one of the two 
quantities in the ratio used to determine n s , extrapolate or interpolate as needed to reach 
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TABLE VI: Values of n s and a -1 (in GeV) for various methods of the analysis that do not require 
(long) chiral extrapolation. In each box, the upper entry is from set R; the lower, from set P. 
These should be compared with our central values, which come from linear chiral extrapolation 
of f qq (and quadratic interpolation of mjL): k s = 0.15873(5), a -1 = 2.19(3) GeV (set R) , and 
k s = 0.16152(14), a -1 = 1.71(4) GeV (set P). We also show the changes (in MeV) from central 
values that each method produces in decay constants, as well as the average and standard devation 



of the mean of those changes. These quantities are used in Sec. VD to estimate the systematic 
effects of chiral logarithms. Lines for which a -1 differs by more than 20% from the central value 
(indicated by '*') are considered unreliable and are eliminated from the averages 



method 


description 




a" 1 


fn 


fn 


fn 


fn 


1 


Kg. TTlggj 'fss : 


■ f ^ ^ 

u " J SS 


15909(7") 


2.38(4) 


+ 11 


+17 


+2 


+6 




777 ■ pn (111 h " 


Jss- c 4- VM-^y 


0.16137(15) 


1.68(4) 


-3 


-7 


+7 


-12 


2 


K, o " 777 I fr.,. ' 
™S • ' 1 fc SS / J SS J 




15926(7") 


2.48(4) 


+13 

\ J.KJ 


+28 


+0 


+13 




779 ■ pn ( LL 1 i • 


f • en fEnh 
Tss • eq. U^JJ 


0.16166(14) 


1.76(4) 


-2 


+3 


+6 


-4 


3 


Kg . TYlss 1 fss ; 


rf7 ' 

u " J SS 


15847(8") 


2.22(3) 


+8 


+5 


+4 


_3 




m • en rt43l"l • 
m ss . cq. \p<J\j , 


f ■ en (1391") 
Jss- e 4- V|23l 


0.16024(22) 


1.55(4) 


-7 


-21 


+7 


-21 


4 


r>>S ■ ln SS/ J SS ) 


Uj - J SS 


15937(7") 


2.47(4) 


+13 


+24 


-Li 

1 -L 


+10 




777, ■ ea. (1121) : 


f QQ ■ en. (B9U 


0.16185(14) 


1.75(4) 


-2 


-1 


+6 


-7 


5 


K s . Tfl ss / f ss , 


a: /ss 


0.15953(7) 


2.57(4) 


+15 


+34 


-1 


+17 




m ss : eq. (gD ; 


/ss : eq. © 


0.16211(13) 


1.83(5) 


-0 


+9 


+5 


+0 


6 


k s : niss/m^ ; 


a: 777^ 


0.15922(6) 


2.46(4) 


+13 


+21 


+1 


+7 




m ss : eq. (|l|) ; 


/ss: eq. (||) 


0.16232(7) 


1.96(2) 


+3 


+11 


+4 


-2 


7 


k s : mss/m^ ; 


a: 7770 


0.15922(6) 


2.46(4) 


+13 


+27 


+1 


+13 




m ss : eq. @) ; 


fss : eq. (0) 


0.16232(7) 


1.96(2) 


+3 


+16 


+4 


+4 


8 


k s : mss/m^ ; 


a: 7770 


0.15871(6) 


2.33(3) 


+10 


+10 


+3 


-0 




m ss : eq. (||) ; 


fss- eq. (H) 


0.16170(7) 


1.89(2) 


+1 


+3 


+5 


-9 


9 


k s : mss/m^ ; 


a: 7770 


0.15946(6) 


2.53(4) 


+14 


+27 


-1 


+11 




m ss : eq. © ; 


fss ■ eq. (H) 


0.16262(7) 


2.00(3) 


+3 


+15 


+3 


+2 


10 


k s : mss/m^ ; 


a: 7770 


0.15946(6) 


2.53(4) 


+14 


+32 


-1 


+17 




m ss : eq. © ; 


/ss : eq. © 


0.16262(7) 


2.00(3) 


+3 


+20 


+3 


+8 


* 11 


K S : fss/ m (f> ! 


a: /ss 


0.16034(107) 


2.83(47) 


+18 


+48 


-5 


+27 


* 


fss ■ eq. 


dD 


0.16479(71) 


2.36(16) 


+9 


+47 


-2 


+35 


12 


K S : fss/ rtl 4> ! 


a: /ss 


0.15881(80) 


2.35(21) 


+11 


+17 


+3 


+6 


* 


/ss : eq. 




0.16420(67) 


2.25(14) 


+7 


+42 


-0 


+31 




average 






+6 


+14 


+3 


+3 




standard deviation of mean 






2 


3 


1 


2 
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FIG. 12: Quadratic extrapolation of m qq / f qq to m ss / f ss on set R, with m ss from eq. ( |4l| ) and /. 



from eq. (39). The abscissa of the burst gives the value of the quark mass, aM^, at k s . This 



particular determination is called "method 1" (see text and Table |Vl| ). 

k = k s , and set the result to the "physical" value of that quantity. The results from either of 



the two quantities in the ratio should be consistent; they are. Figure [14| shows a quadratic 
extrapolation of f qq to k s ; a~ l is fixed from the result via eq. ([39]). Although an extrapolation 
is again required in this particular case, it is only over a short distance in quark mass. 

Once the scale is determined, the standard extrapolation of rr? produces the light quark 
hopping parameter, k u 4. It is, of course, very close to n c in all cases. 



Results for various scale determinations are shown in Table [VI]. In most cases, the val- 
ues of a -1 there are significantly larger than those from the standard linear (or quadratic) 
extrapolation of f n (see Table [V]). This is not unexpected because extrapolation from rela- 
tively large mass without the chiral log term in eq. should overestimate af n . Further, 
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FIG. 13: Quadratic interpolation of rn qq /m^ to m ss /m^ on set R, with m ss from eq. The 
abscissa of the burst gives the value of the quark mass, aM 9i 2, at k s . This particular determination 



is called "method 8" (see text and Table VI). 



other lattice spacing determinations from light quark physics (e.g., m p , mjv) also typically 
involve linear or quadratic chiral extrapolations from rather high masses, so their generally 
good agreement with the standard f n scale cannot be used to rule out the results in Ta- 
ble [V[ Indeed, it has been known for some time that heavy quark physics (charmonium or 
upsilonium) typically gives scales ~ 10 to 20% larger than light quark physics. (For a recent 
example, see Table I in Ref. f5lf .) It now appears that at least some of this discrepancy is 
due to the extrapolation in the light quark mass. 

On the other hand, the result for a" 1 is unreasonably large (>30% bigger than the central 
value) in three cases in Table [VT|. (In all other cases the scale is at most 17% greater than 
the central value.) The three cases also have long extrapolations to find value of 
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FIG. 14: Quadratic extrapolation on set R of f qq to n s , which in turn was found by extrapolation 



of m ss / f ss using eqs. (|42| ) and (|38|). This is method 4 in Table [Vl 



k s very close to n c — compare Table |III|) and very large statistical errors in both k s and a . 



These cases are marked with asterisks and are omitted from any averages of decay constant 
effects. 

Given k s , n U) d, and a" 1 , the ratio f gq / fq g can now be interpolated/extrapolated as in 
Fig. [XT] to k s and k u ^- Using the physical f n and the relevant choice for / ss , this produces 
/q s or fQ U ^di which are then interpolated in heavy quark mass in the same way as in our 
standard analysis, described in the next section. The differences of the final decay constants 



from the central values are displayed in Table VI for each of the methods. These differences 



will be used in Sec. [VD| to estimate the effects of chiral logs in the dynamical case. We also 
explain there why we think it would be inappropriate at this stage to correct the central 
values by the chiral log effects. Instead we use the changes shown in Table [Vj only to 
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estimate the systematic error. 



F. Interpolations in Heavy Quark Mass 

We proceed to compute physical decay constants such as /# and fs s for each lattice set. 
Our starting point are the values of fq q Mq q (for q = u, d or q = s) as functions of the 
heavy-quark mass produced by the chiral fits of the previous two sections. 

The static limit is also included where we have it. According to the heavy quark effective 
theory (HQET) [j5^| , fq q «/ Mq 9 should depend on the heavy meson mass as a polynomial in 
1/Mq q , up to logarithms. We therefore first divide out the one-loop logarithmic dependence 
of the decay constants in the heavy quark limit |53fl , producing what we call f'Q q \J Mq q : 



= l + a!{q^ q M Qq )/, ' (44) 

where we have ignored the difference between the heavy quark and heavy meson masses, 
and where q* takes the values discussed in Sees. [11 A| and 11 B . 



The data is now expressed in physical units, always using /„. to set the scale for the 
central values. The quantity /A Mq q is then plotted vs. l/Mgg, where Mq q is the kinetic 
meson mass MQ g>2 defined in eq. (|§|). We fit to a polynomial in 1/Mq q , interpolate to m B , 
m B a , m D or m£) s , and then replace the logarithm in eq. ([44]), evaluated at the appropriate 
meson mass. These are always interpolations, not extrapolations, because we have either 
the static-light point (all Wilson sets) or heavy-light masses above the B (the clover sets 
CP1 and J, using the Fermilab formalism). 

For the quenched Wilson data, we do two versions of the polynomial fit: (1) a quadratic fit 
to heavy-light mesons in the approximate mass range 2 to 4 GeV plus the static-light meson 
( "heavier- heavies" ) and (2) a quadratic fit to mesons in the approximate mass range 1.25 to 
2 GeV plus the static-light meson ( "lighter- heavies" ) . These fits keep just one eigenvector 
of the correlation matrix, which corresponds to A cut = 0.9 to 1.1. To estimate the effect 
of leaving out higher powers in 1/Mq q in the fits, we also perform, for the central q* and 
chiral fit choices, fit (3): a cubic fit to all the mesons in the range 1.25 to 4 GeV plus the 
static-light meson. The correlation matrix for fit (3) typically has almost twice the number 
of eigenvectors as fits (1) and (2), and we keep 2 of the them. This corresponds to A cut = 0.2 
to 1.0. 
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We make basically the same fits for the Wilson data on the dynamical lattices. The main 
difference is that we cut off fits (1) and (3) at approximately 3, rather than 4, GeV. These 
lattices are almost all quite large, and, as explained in Sec. [Ill A| , we have trouble pulling out 
the lightest state for very heavy masses on large lattices with our approach to the hopping 
expansion. To make up for some of the points lost by the reduced upper cutoff on fit (1), 
we also reduce the lower cutoff slightly, to 1.8 GeV. 

For the quenched clover sets, we make corresponding fits. However, the mass ranges are 
somewhat different because we have only five heavy quark values, do not have a static point, 
and, most importantly, use standard algorithms with FFT, facilitating the extraction of the 
lowest states even for very heavy masses. In this case, fit (1) (heavier-heavies) is a quadratic 
fit over the approximate meson mass range 2.3 to 6 GeV; while fit (2) (lighter-heavies) is 
over the mass range 1.7 to 3 GeV. 

In central values, we use fit (1) for and fs s and fit (2) for fo and fjj 3 . The alternative 
fits go into the systematic error estimates, as in Ref. [|J. However, for the central values of 
ratios involving both B and D mesons (i.e., fs/ /d s , /b s //d s , and Jb/ Jd), both numerators 
and denominators are taken from fit (2). As explained in Sec. [TV], this tends to reduce the 
estimate of the magnetic mass error. 

Figures |TS], |TB|, and [H] give examples of the behavior of fQ q ^Mg q for the quenched 
Wilson, Wilson on dynamical lattices, and quenched clover cases, respectively. Essentially 
all fits on all sets are acceptable. In Fig. |TB|, one can see the tail-off of fQ q \/MQ g for large 
heavy quark masses (0.2 GeV -1 < 1/Mq q < 0.3 GeV -1 ). As mentioned above, we attribute 
this to contamination by higher momentum states, which, for large masses and volumes, are 
very close in energy to the zero momentum state. These points are therefore not included 
in the fits. Note that the term of order 1/Mq is not reliably determined in our data; it 
changes sign between the Wilson and the clover cases. This is not surprising since in neither 
case is the formalism correct through order 1/Mq q for aMg ~ 1. 

Using fits like those in Figs. pT5HT7|, we now interpolate the data, replace the perturbative 
logarithm in eq. (|44"|), and divide by the appropriate 'Mq q to find fs, fB s , Jd, and f Ds 
for each data set. The resulting decay constants and ratios will be extrapolated to the 
continuum in Sees [TV] and [V]. Before doing so, however, we repeat the analysis so far for all 
the other ~ 25 versions of reasonable plateau choices, as discussed in Sec. |111 Q We then 
find the standard deviation of the results over the other versions and add it in quadrature 
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FIG. 15: f Qq \jM Qq vs. l/M Qg for set CP (quenched Wilson). The scale is set by f n . The solid 
line is fit (1) ("heavier-heavies") and includes points marked with a cross. The dotted line is fit 
(2) ("lighter-heavies") and includes points marked with a plus. The fits have Xcut/d-o.f. = 0.4 and 
0.9, respectively, with A cut = 1 (1 eigenvector kept). 

with the raw jackknife error of the central value. Henceforth, the statistical error of any 
quantity will be taken from the result of this procedure. Typically the procedure increases 
the statistical errors by ~ \^2; we believe it mitigates any biases introduced from our choice, 
for the central values, of the fits with lowest statistical errors and x^ ut /d.o.f . 



In Table [VI]] , we collect the central values of fs, /b 3 , Id, and fn s for the various sets. 
Similarly, Table |Vlli| gives central values of ratios f Ba /fB, /d s //d, /b//d s , /b s //d s , and 
Ib/ Id- 
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FIG. 16: Same as Fig. 15, but for set R (Wilson valence quarks on dynamical lattices). The solid 
line has x^ ut /d.o.f. = 0.2; the dotted line, 0.3. The tail-off of fn q \/MQq in the range 0.2 GeV- 1 < 
1/Mq q < 0.3 GeV -1 is attributed to the difficulty in isolating asymptotic states for large masses 
and volumes — see text. 

IV. QUENCHED APPROXIMATION RESULTS 

Final results and errors in the quenched approximation are determined much as in Ref . [|J . 
However, there are some significant differences, especially for the continuum extrapolation 
and the estimate of the associated errors. We discuss our methods in detail where they 
differ from where the methods are the same, we include only a very brief description for 
completeness. 

We begin with the continuum extrapolation of various quantities. We focus on f B , f Bs , 
fB a l fs, and fo s , which are probably the most important, phenomenologically. Figures |18| - 
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FIG. 17: Same as Fig. 15, but for set J (quenched clover), with NP-tad renormalization. The solid 



line has Xp Ut /d.o.f. = 0.55; the dotted line has no degrees of freedom. The solid line is slightly 



concave down, unlike the case in Figs. 15 and 16 



21| show the data for these quantities as a function of lattice spacing. The behavior of the 
other decay constants and ratios is similar. 

It is not a priori obvious how to extrapolate decay constants and ratios to the continuum. 
As discussed in Sec. [11 A| , our Wilson valence results have a v errors as well as errors of the 
form ah.QCD x /i(oMq). Here h(aMo) is a calculable (in perturbation theory) function that 
is expected to be 0(1) everywhere. 9 Since aMq>l for our entire range of a values, the 
assumption of a dominantly linear dependence on a is only one possibility. A practical 



9 Despite the fact that static quarks are trivially 0{a) improved, the function h(aMq) does not vanish even 
as Mq — > oo for fixed a, because the Wilson light quarks still have 0(a) errors. 
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TABLE VII: Central values of decay constants, in MeV, for each data set. Statistical errors include 
the effect of changing the time ranges over which correlators are fit, as described in the text. For 
sets A and B, the values reported are those for which the light-light results (k c , k s and a -1 ) and 
their errors are taken from the averages over set 5.7-large. 



name 


fB 


fs s 


fD 


f'D s 


quenched Wilson 


A 


193.6(8.6) 


234.9(5.8) 


216.6(7.8) 


256.9(6.5) 


B 


196.8(11.1) 


236.1(8.2) 


220.7(10.9) 


261.8(7.4) 


E 


190.9(7.7) 


219.1(7.8) 


214.7(8.6) 


246.1(8.1) 


C 


172.3(8.0) 


206.2(6.8) 


198.8(7.7) 


232.8(6.9) 


CP 


177.0(7.8) 


210.8(6.4) 


206.7(6.4) 


238.2(5.7) 


D 


174.9(7.5) 


199.8(6.5) 


206.8(7.9) 


232.5(6.3) 


H 


180.6(12.1) 


206.7(10.8) 


206.6(10.9) 


232.9(9.0) 




N f 


= 2 Wilson 






L 


188.6(9.7) 


220.3(10.3) 


214.7(7.4) 


249.5(6.7) 


N 


205.7(13.6) 


239.0(10.4) 


222.9(10.6) 


261.5(7.6) 





206.8(12.6) 


239.9(10.8) 


230.8(6.3) 


262.4(5.0) 


M 


190.6(12.5) 


226.9(10.2) 


215.9(11.5) 


250.2(8.4) 


P 


193.1(6.9) 


225.6(6.4) 


212.9(6.7) 


249.5(5.9) 


U 


196.5(9.4) 


235.0(7.9) 


224.8(7.2) 


261.2(7.6) 


T 


193.3(15.8) 


219.3(12.5) 


209.1(8.3) 


236.3(6.7) 


S 


202.6(6.8) 


234.6(5.6) 


223.9(5.2) 


256.6(4.7) 


G 


198.5(6.2) 


234.0(6.7) 


220.0(5.0) 


254.7(5.2) 


R 


206.2(7.9) 


239.2(8.1) 


223.4(5.4) 


254.6(5.1) 


quenched clover 


CP1.NP-IOY 


184.1(5.7) 


212.8(4.4) 






CPl_NP-tad 


176.1(5.2) 


203.4(3.9) 


196.3(3.7) 


220.0(2.8) 


J.NP-IOY 


176.6(6.3) 


204.4(5.9) 






J_NP-tad 


174.0(6.0) 


201.9(5.7) 


203.5(4.8) 


228.2(4.1) 



alternative is the assumption that, for a smaller than some value, the errors are small 
enough that the difference with continuum values is negligible — so that extrapolation with 
a constant function is warranted. 



We confront these assumptions with the data in Figs. 18-21. For the Wilson valence 
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TABLE VIII: Same as Table 



VII 



, but for ratios of decays constants. 



name 


fB a 1 !b 


fDs 1 /d 


fB/fD 3 


/b s //d s 


Ib/Id 


quenched Wilson 


A 


1.213(39) 


1.186(21) 


0.768(24) 


0.926(20) 


0.911(27) 


B 


1.200(30) 


1.187(29) 


0.793(31) 


0.928(28) 


0.941(50) 


E 


1.147(34) 


1.146(16) 


0.739(40) 


0.916(33) 


0.846(46) 


C 


1.197(22) 


1.171(17) 


0.740(28) 


0.888(25) 


0.867(32) 


CP 


1.191(21) 


1.152(12) 


0.752(36) 


0.899(26) 


0.867(41) 


D 


1.142(16) 


1.124(16) 


0.755(23) 


0.866(20) 


0.849(26) 


H 


1.145(22) 


1.128(18) 


0.776(28) 


0.888(19) 


0.875(25) 




N f ~- 


= 2 Wilson 






L 


1.168(16) 


1.162(16) 


0.777(30) 


0.918(22) 


0.903(38) 


N 


1.162(33) 


1.173(25) 


0.788(37) 


0.927(29) 


0.924(41) 





1.160(25) 


1.137(18) 


0.789(32) 


0.916(23) 


0.897(30) 


M 


1.191(30) 


1.159(26) 


0.793(25) 


0.928(22) 


0.920(25) 


P 


1.168(12) 


1.172(12) 


0.784(16) 


0.914(12) 


0.919(18) 


U 


1.196(28) 


1.162(12) 


0.763(41) 


0.904(37) 


0.886(42) 


T 


1.134(34) 


1.131(16) 


0.840(41) 


0.958(32) 


0.950(41) 


S 


1.158(17) 


1.146(11) 


0.790(19) 


0.927(15) 


0.905(19) 


G 


1.179(14) 


1.158( 8) 


0.788(20) 


0.937(26) 


0.912(22) 


R 


1.160(15) 


1.140( 8) 


0.814(28) 


0.947(20) 


0.927(32) 



quenched clover 



CP1_NP-I0Y 


1.156(16) 


1.120( 8) 








CPl_NP-tad 


1.155(16) 


1.121( 9) 


0.800(19) 


0.923(12) 


0.896(18) 


J.NP-IOY 


1.157(14) 


1.117( 8) 








J_NP-tad 


1.160(14) 


1.121( 8) 


0.772(21) 


0.892(18) 


0.866(24) 



data, we show linear fits over all a and constant fits for a < 0.5 GeV~ (/3 > 6.0). Both 
types of fits are generally quite good. The exception is the constant fit for the ratio f B J f B . 
In H, the relatively poor confidence level of the constant fit for f Bs /f B (or relative 
to that of the linear fit led us to choose a linear extrapolation for the central value of the 
ratios. That, in turn, required choosing the linear extrapolation for the central values of 
the decay constants themselves, since it would be inconsistent to assume linear behavior 
for f B J ' f B but constant behavior for f B and f Ba separately. Note, however, that if we just 
look at the two finest lattices (a < 0.36 (GeV) -1 , (3 > 6.3), the behavior of f B J f B is quite 
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FIG. 18: /b vs. a for quenched lattices; the scale is set by /„■. Diamonds are results with Wilson light 
quarks and Wilson or static heavy quarks. Octagons and crosses are results with nonperturbative 
clover heavy and light quarks; "NP-IOY" (octagons) and "NP-tad" (crosses) differ in how the 
renormalization of the heavy quarks is performed (see text). For clarity, the octagons have been 
moved slightly to the right, and the fit to the crosses has been slightly lowered. 



consistent with a constant; such a fit is also shown in Fig. 



The new quenched clover data, shown in Figs. |18|-|21| for both the NP-IOY and NP-tad 
schemes, has clarified the situation somewhat. The discretization errors here should be 
considerably smaller than for Wilson valence quarks. As discussed in Sec. [II B| the errors are 
formally 0(a 2 A 2 QCD ) and either 0{a 2 v ) (NP-IOY), or Q{a 2 M 2 Q ) (NP-tad). Because there will 
be also be a function like /i(aMg) in this case, the actual behavior with a when oMq ~ 1 
is likely to be complicated. The best we can do with just two clover data points is to 
assume that the errors are small enough that a constant extrapolation is warranted; such 
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FIG. 19: Same as Fig. [D], but for fs s - For clarity, the octagons have been moved slightly to the 
right. 



fits are shown in Figs. |l8|-pl|. Comparable extrapolation of clover data with a constant was 



performed in Refs. [17, 54 



For f B J f B the clover data show very little a dependence and give a result compatible with 
the various constant fits to the small-a Wilson data. The clover results are not compatible 
with the linear extrapolation of the Wilson data, which is now seen to give a rather low result. 



Recent preliminary quenched results [45| with clover valence quarks on Symanzik improved 
glue are also incompatible with the Wilson linear extrapolation. For our central quenched 
value of fsj Jb or foj 1 fo we therefore drop the linear Wilson extrapolation and average the 
four constant extrapolations: two for Wilson (a < 0.5 GeV -1 and a < 0.36 GeV -1 ) and two 
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FIG. 20: Jb 3 / Ib vs. a for quenched lattices. Labels are the same as in Fig. [TJ|, but one additional 
fit is shown: a constant fit to the two diamonds (Wilson quark results) with smallest a (a < 
0.36GeV _1 ). For clarity, the octagons have been moved slightly to the right and the fit to the 
octagons has been slightly lowered. 

for clover (NP-IOY and NP-tad). 10 The systematic error of the continuum extrapolation is 
then taken as the standard deviation of the four individual extrapolations. The other decay 
constant ratios (/b//d s , /b s //d s , and /b//d) are treated similarly, although there is one 
fewer result to average, since NP-IOY is not applicable. 

Though we have dropped the linear extrapolation from the analysis of the ratios, it is not 

10 Although the calculation in Ref. [^8| is not well controlled at the D mass, the NP-IOY procedure may be 
used for JdJ Id because the renormalizations cancel. Note that NP-tad involves the light quark mass in 
the renormalizations (see eqs. |l^ and |l9|) so does not give identical results to NP-IOY even for JbJ Jb or 
IdJ Id- 
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FIG. 21: fn s vs. a for quenched lattices. Labels are the same as in Fig. 18, but the NP-IOY points 
have been omitted because the perturbative calculation is not available at the relevant lattice 
masses. 



inconsistent to include it in the analysis of the decay constants themselves. Indeed, for fg 
and Jb s , the downward trend of the clover data as a decreases makes it difficult to rule out 
the linear extrapolation of the Wilson data. On the other hand, constant clover extrapola- 
tions do give results closer to the constant Wilson extrapolations than to the linear Wilson 
extrapolations. For fjj a (and fo, not shown), the situation is reversed: The clover data has 
an upward trend as a decreases; yet constant clover extrapolations give results (slightly) 
closer to the linear Wilson extrapolations than to the constant Wilson extrapolations. To 
obtain the central values of the decay constants, we therefore average the results of all the 
extrapolations and take the standard deviation of the results as the continuum extrapola- 
tion error. For fg and fs 3 a total of four fits are included: linear Wilson, constant Wilson, 
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constant NP-IOY and constant NP-tad. For fn and fo s there are three fits, since NP-IOY 
is omitted. 

As described in Sec. Q I fD| , we estimate the chiral extrapolation errors by comparing 
(after continuum extrapolation) the central values (which use "chiral choice I") with those 
obtained by changing the chiral fits of the heavy-light and light-light decay constants from 
linear to quadratic ("chiral choice II"). 11 The errors of the chiral extrapolation and other 
systematic errors within the quenched approximation are collected in Tables [IX| and |X|. 
Note that the quoted chiral errors are all positive. This can be traced to the effect of the 
quadratic extrapolation of f qq (used to set the scale through /„.), which is clearly, though 
slightly, concave down (see Figs. || and |7j). The concavity in fq q in the region of the B is 
less pronounced. 



TABLE IX: Central values (f n scale) and errors in MeV for the quenched decay constants. The 



statistical errors and the effects of excited states are combined, as described at the end of Sec. 1IIF . 
Errors marked with explicit + or — signs are treated as signed; all others are treated as symmetric. 
The scale and k s errors are not included in the total error within the quenched approximation but 
are shown for completeness. 





fs 


Ib s 


fn 


Id s 


central value 


173.0 


198.8 


199.5 


223.2 


errors 


statistics & excited states 


5.7 


4.7 


5.6 


4.6 


continuum extrapolation 


8.7 


14.8 


4.9 


11.1 


chiral extrapolation 


+8.6 


+9.4 


+4.3 


+6.5 


perturbative 


9.6 


14.1 


6.6 


10.9 


magnetic mass 


1.7 


1.9 


5.1 


5.7 


1/M fit 


2.0 


1.7 


0.5 


0.2 


finite volume 


+2.8 -8.5 


+1.0 -7.8 


+3.8 -5.1 


+4.3 -4.1 


scale (change to m p ) 


-3.5 


-5.9 


+4.2 


+4.0 


k s (change to <p) 




+3.7 




+2.3 



In the dynamical case, we attempt to estimate an additional chiral error coming from chiral logarithms 
by performing a separate chiral extrapolation of f qq / fcjq- (See Sec. HIE and VD.) This is not feasible in 
the quenched case: the quenched chiral logs in /q 9 have coefficients with unknown magnitude and sign 
Hi ; while f qq has no quenched logs at all at one loop Jsg) . 
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TABLE X: Same as Table IX but for decay constant ratios. 





Jb s 1 Ib 


Id s / Id 




fB s /fD s 




central value 


1.155 


1.128 


0.769 


0.891 


0.871 


errors 


statistics & excited states 


0.011 


0.008 


0.015 


0.012 


0.016 


continuum extrapolation 


0.009 


0.012 


0.017 


0.019 


0.013 


chiral extrapolation 


+0.003 


+0.014 


+0.009 


+0.009 


+0.020 


perturbative 


0.008 


0.011 


0.015 


0.020 


0.017 


magnetic mass 


0.000 


0.000 


0.015 


0.018 


0.017 


1/M fit 


0.001 


0.000 


0.006 


0.008 


0.009 


finite volume 


+0.012 -0.013 


+0.008 -0.000 


+0.025 -0.000 


+0.003 -0.009 


+0.028 -0.000 


scale (change to m p ) 


+0.001 


-0.001 


+0.009 


+0.008 


+0.007 


k s (change to (f) 


+0.025 


+0.018 


-0.009 


+0.004 





The perturbative error is estimated by varying, over a "reasonable range," the values of 
q* used in the one-loop renormalization constants. For Wilson fermions, we take the range 
for the heavy-light currents to be 1/a < q* < 2.86/a, with 1.43/a the central value, as 
described in Sec. [LI Aj . Similarly, for the light-light Wilson currents q* ranges between 1/a 
to 4.63/a, with 2.32/ a the central value. In the clover case, perturbation theory for the 
heavy-light currents is only relevant for NP-IOY. For central values, we take q* = 3.34/a 
(set CP1) and q* = 2.85 /a (set J), which come from the static-light calculation of [|TT[] with 
the corresponding clover coefficients. The scale q* is then allowed to range between 1/a and 
twice the central value. For light-light clover currents, only b A is treated perturbatively; the 
central value for q* is taken to be 1/a (see Sec. |IIB| ). This gives the central values for b A 
shown in Table |1[ The upper end of the range of b A shown comes from taking q* = 0.7/a; 
the lower end, from using "boosted perturbation theory" with g 2 = 6/(j3 < P >). (< P > 
is the mean plaquette, normalized to have maximum 1). This is equivalent to taking a q* of 
roughly 5.25/ a, so we are using a rather conservative range. 

As mentioned in Sec. |II A , there is a systematic error associated with the fact that 



c mag = M 2 /M 3 is not equal to 1 with Wilson fermions. Because c mag has complicated 
dependence on a, this error is not removed by any of the simple extrapolations available to 
us. One may argue that the residual effect is just one particular discretization error and 
therefore has already been included. However, if one models this error for both linear and 
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constant extrapolations using eq. (f|) for Mi and M3 (along the lines of what was done in 
Refs. JIT) and Q), one finds that the error is larger with a constant extrapolation but has 
the same (but unknown) sign in both cases. Therefore we believe it reasonable to include 
as an additional error the linear extrapolation estimate of the Wilson magnetic mass error. 
From the tadpole improved tree-level model, one estimates these errors as ~2% for fs and 
~3% for /d (see ||). An alternative estimate comes from the comparison of the results of 
interpolations to the physical heavy meson masses using the "heavier-heavies" (fit (1) - see 
Sec. 1 1 IF) with those using the "lighter-heavies" (fit (2)): the lighter masses are affected 
much less by the magnetic mass error, and the static point is not affected at all. We take 
the larger of the two estimates as our magnetic mass error for Wilson fermions. 

The magnetic mass error is absent for clover fermions. Therefore, in our final error budgets 
(Tables [IX|, |X|) we multiply the Wilson magnetic mass error by 1/2 or 2/3, depending on 
the relative number of Wilson and clover estimates that go into the central value. 

Note that the magnetic mass errors in the tables are considerably smaller for B mesons 
than for D mesons, despite the fact that the difference between M2 and M3 increases with 
the lattice mass. The point is that the magnetic mass errors are systematic effects on the 
1/Mq corrections, and such corrections are inherently bigger for D mesons than for 5's. 
Further, especially large errors can be introduced if 1/Mq fits in the range of the B are 
extrapolated back to the D region. For that reason we always use fit (2) ("lighter-heavies") 
for central values of ratios that involve both D's and B's: fs/ fD s , fB s / fo s , and Jb/ Id- 

The remaining two systematic errors, the effect of the interpolation in l/MQ q and the 
finite size errors, are estimated just as in Refs. |Q. For the central values, we truncate the 
fit of fQ q a/ Mq q vs. 1/Mq q at quadratic order. We estimate the error thereby introduced by 
changing to cubic fits (with mass range 1.25 to 4 GeV, plus the static point when available). 
The errors found are ~ 1 %; this is what one would expect if the mass scale of the cubic term 
is ~0.75 GeV, roughly the scale size found in the linear and quadratic terms. 

We estimate the finite volume effects by finding the fractional difference between results 
on set A (spatial size ~ 1.2 fm) and set B (~ 2.5 fm). Since set A is smaller than the 
other quenched lattices (~ 1.3-1.5 fm) and B is much larger, this should bound the finite 
volume error. To be conservative, we consider both: (a) the difference when all quantities are 
computed individually on sets A and B and (b) the difference when the light-light quantities 
are held fixed to their values from set 5.7-large. Since /„- generally suffers larger finite size 
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effects than fg q , these two estimates typically have opposite signs; in that case we include 
both estimates as signed errors. When the estimates have the same sign, however, we simply 
choose the larger. 

Tables [TX| and |X] also show errors associated with fixing the scale (changing from /„. to 



nip) and fixing k s (changing from using the pseudoscalars to using the 4> meson). Logically, 
these should be considered errors of the quenched approximation, not within the quenched 
approximation, and are not included in this section. Indeed, the question "what is /# in 
the quenched approximation?" is only well defined when one specifies how the scale is fixed. 
Even in the continuum limit, different scale choices (and different ways of fixing k s for 
strange-quark quantities) must give different results in the quenched approximation. The 
differences should of course go away in the continuum limit of the full theory. In Sec. |V], 
where we attempt to quote results that can be directly compared with experiment, such 
errors are taken into account. 

Our final results for heavy-light decay constants within the quenched approximation (fix- 
ing the scale by /„.) are 

f B = 173(6)(16) MeV; f Bs = 199(5)(^) MeV 
Id = 200(6)(11) MeV; f Ds = 223(5)(±#) MeV 
^ = 1.16(1)(2); ^ = 1.13(1)(2) (45) 

JB JD 



p- = 0.77(2)(±£); = 0.89(1 w+4 

JDs ' jD a 



^ = 0.87(2)(H) . 

JD 

The errors are statistical and systematic (within the quenched approximation), respectively. 
Relevant systematic errors in Tables [IX| and |X| have been combined in quadrature. Errors 
whose signs are not likely to be reliably determined by our procedures (continuum extrapo- 
lation, perturbation theory, magnetic mass, 1/M fits) have been treated as symmetric errors. 
The others (chiral extrapolation and finite volume) have been treated as signed errors. The 
results in eq. (|45"D differ from those in Ref. Q due to: 



1. Inclusion of new data from sets CP, CP1 and J. 

2. Setting the central value of the heavy-light scale from the static-light calculation 
of Ref. 0, rather than that of Ref. J2|. 
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3. Other changes in analysis, motivated by the new runs. The most important of 
these is the way we find the central value of the continuum extrapolation (as discussed 
above, we now average our four possible versions rather than taking only the linear 
Wilson fit). In addition, the details of the error estimate for the chiral extrapolation 
have changed. Some alternative chiral fits used previously — e.g., linear fits of vs. 
quark mass — are convincingly excluded by the new data. 



V. RESULTS WITH DYNAMICAL QUARKS 
A. Continuum Extrapolation 

Dynamical Nf = 2 results for fs, /r,, JbJ Ib, and fu a as a function of lattice spacing are 



shown in Figs. |22], |23|, and respectively. Leaving aside the "fat clover" results for now, 
the data in all cases seem to favor constant fits; indeed, the best linear fits have very small 
slopes. Note however that the smallest lattice spacing here is ~ 0.45 (GeV) -1 ~ 0.09 fm; 
whereas in the quenched case we have data down to ~0.23 (GeV) -1 ~ 0.045 fm. It is thus 
possible that the apparent independence of lattice spacing is due to the cancellation of two 
effects: (1) an overall decrease as lattice spacing decreases, which was one of the alternatives 
considered in the quenched case, and (2) the turning on of short distance dynamical fermion 
effects as one moves away from the quite coarse spacings of sets L and N. The latter effect 
could be exacerbated by staggered flavor violations, which would be especially large on the 
coarsest lattices and which would reduce the effective number of dynamical flavors. 

If the above possibility is realized, then the Nf = 2 results could well begin to decrease for 
still smaller lattice spacings as the quenched-like behavior sets in. For the decay constants, 
we therefore consider two alternative extrapolations: the constant extrapolation of all Nf = 2 
data, and a linear extrapolation that begins at the average value of the results on the two 
finest lattices (sets R and G) and then continues to the continuum limit with the quenched 
slope (see Figs |22|, [23|, and ES[) For ratios of decay constants, we ruled out the linear 
extrapolation in the quenched case. Yet the two finest quenched lattices (D and H) have in 
general lower values for the ratios than the averages that include the quenched sets (C and 
CP) that are comparable to the finest Nf = 2 lattices. The two alternatives for ratios are 
therefore taken to be (1) the constant extrapolation of all Nf = 2 data, and (2) the first 
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extrapolation reduced by the quenched difference: (average of C, CP, D, and H) - (average 
of D and H). Figure p4] shows these alternatives. In all cases we then take the central value 
to be the average of the two alternatives, and the error of the continuum extrapolation to 
be the "sample standard deviation" of the two (dividing by n — 1 = 1, not n — 2). Central 
values and errors for the Nf = 2 data are shown in Tables |X| and |XTT. 



TABLE XI: Central values (f n scale) and errors in MeV for the dynamical (Nf = 2) decay constants. 
As in Tables IX and [X|, the statistical errors and the effects of excited states are combined. The 
errors above the line (i.e., up to and including finite volume errors) are treated as errors within 
the Nf = 2 partially quenched approximation. Errors below the line are treated as errors of that 
approximation. In general, errors marked with explicit + or — signs are treated as signed, and 
other errors are treated as symmetric. The exception is partial quenching, where we do not take 
the sign seriously but show it nevertheless in parentheses. 





f B 


fBs 


Id 




central value 


190.5 


217.3 


214.9 


241.0 


errors 


statistics & excited states 


7.1 


6.4 


6.1 


5.2 


continuum extrapolation 


11.3 


21.0 


8.5 


18.7 


valence chiral extrapolation 


+16.6 


+14.7 


+7.5 


+8.3 


perturbative 


12.0 


18.5 


8.2 


15.1 


magnetic mass 


3.8 


4.4 


8.4 


9.7 


1/M fit 


2.6 


2.7 


1.0 


0.9 


finite volume 


+7.7 0.0 


+5.2 -> 0.0 


+3.4 -> 0.0 


-0.1 0.0 


partial quenching 


(+)2.4 


R3-0 


(+)3.4 


R3-8 


scale (change to m p ) 


+10.6 


+8.7 


+5.4 


+3.9 


k s (change to <j>) 




+3.9 




+2.3 


missing dynamical s quark 


+8.7 


+9.2 


+7.7 


+8.9 



B. Fat Links 

In the above discussion of the continuum extrapolation, we ignored the fat-link clover 
results. If taken at face value, these results would imply the existence of extremely large 
discretization errors. We therefore need to examine the fatdink computations in detail. 
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TABLE XII: Same as Table XI but for decay constant ratios. 





/b s / Jb 


fD s 1 /d 


fB/fD s 




Ib/Sd 


central value 


1.158 


1.142 


0.793 


0.922 


0.913 


errors 


statistics & excited states 


0.011 


0.009 


0.016 


0.013 


0.016 


continuum extrapolation 


0.015 


0.014 


0.005 


0.004 


0.001 


valence chiral extrapolation 


-0.016 


+0.005 


+0.032 


+0.019 


+0.037 


perturbative 


0.012 


0.011 


0.034 


0.043 


0.042 


magnetic mass 


0.003 


0.002 


0.024 


0.028 


0.027 


1/M fit 


0.001 


0.000 


0.014 


0.020 


0.015 


finite volume 


-0.019 -> 0.000 


-0.018 


+0.026 -> 0.000 


+0.010 -> 0.000 


+0.016 -» 0.000 


partial quenching 


(-)0.023 


(-)0.026 


(+)0.027 


(-)0.021 


(+)0.008 


scale (change to m p ) 


-0.010 


-0.005 


+0.015 


+0.017 


+0.015 


k s (change to <p) 


+0.014 


+0.017 


-0.008 


+0.004 




missing dynamical s quark 


+0.001 


+0.007 


+0.012 


+0.015 


+0.021 



These computations use valence quarks — both heavy and light — with the standard clover 
action, but with gauge links that have first been "fattened" by iV iterations of APE smearing 
The coefficient of the sum of the staples is c/6 and that of the forward link is 1 — c; 
a projection back into SU(3) is included after each smearing step. The fat-link results (set 
RF) displayed in Figs. |2^|-f2o] have N = 10, c = 0.45. The clover coefficient csw is taken 
to have its tree-level value (1.0); this is also approximately the tadpole-improved value, 
since the fattening strongly suppresses tadpole contributions. Physically, APE smearing 
corresponds roughly to a Gaussian smearing of the fermion-gauge field interaction over a 
range (x 2 ) ~ cN/6 ||. 

Various kinds of fat links have come to play a major role in lattice simulations in the last 
few years. The motivation for introducing them in the context of Wilson-like fermions f59fl 
was that they improve the chiral properties of the fermions. This happens in several (related) 
ways: First, fat links reduce additive mass renormalization. They also suppress exceptional 
configurations, which present a severe challenge to clover computations on our dynamical 
lattices |12|, |6(| . (This occurs because they shrink the range of the real eigenmodes of the 
Dirac operator.) Finally, in perturbation theory, fat links bring the vector and axial vector 
renormalization constants Zy and (as well as the scalar and pseudoscalar renormalization 
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FIG. 22: fs vs. a for dynamical Nf = 2 lattices; a few points have been moved a slight distance 
horizontally for clarity. Squares are results with light Wilson valence quarks and Wilson or static 
heavy valence quarks. From left to right, the squares come from (3 = 5.6 (sets G, R, S, T, U), 
(3 = 5.5 (sets P, M, O, N), and (3 = 5.445 (set L). The solid line is a fit of all the Wilson results 
to a constant. The dashed line shows what would happen if the dynamical results decreased for 
smaller lattice spacing with the same slope as the linear fit to the corresponding quenched data. 
The fancy plus is the result with fat-link clover valence quarks (light and heavy) on set RF. The 
fancy cross shows the "corrected" value (see text). The fat clover data (corrected or uncorrected) 
is not included in our final results. 

constants) closer together. 

Simulations of light quark systems with a variety of fat link actions at lattice spacings 
in the range 0.1-0.2 fm show little dependence of physical observables on the amount of 
fattening, even for the very aggressive amount of fattening of the simulations we report 
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FIG. 23: Same as Fig. 22, but for f B 



here. For many quantities, this amount of fattening also gives quite small discretization 
errors 
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We take the light-light renormalization coefficients for fat-link clover fermions from the 



perturbative calculations of Ref . |Tl]] . The heavy-lights (for which perturbative calculations 
do not exist) are normalized using the static-light results of ||11|| . Although, one expects that 
this should be roughly correct for the large values of aM at the B meson, it introduces a 
possibly serious source of systematic error into the fat-link results. 

As first reported in Ref. |T2| , the fat-link clover results for decay constants are seen to be 



much smaller than the apparent continuum-limit results of the Wilson quarks. Simulations 
of QQ systems with fat link quarks also show that fattening suppresses the magnitude of 
vector-pseudoscalar mass splitting. A measurement of the heavy quark potential gives some 
qualitative understanding of both effects: the attractive short distance piece of the potential 
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FIG. 24: /b s / Jb vs. a for dynamical Nf = 2 lattices. Labels are the same as in Fig. 22, but no 
correction to the fatlink clover result is needed for the ratio of decay constants. The alternative 
dashed line assumes a drop when a —* that is the same as the difference in the quenched case 
between a constant fit to the results from the highest three (3 values and a constant fit to those 
from the highest two (3 values. 

is washed away by the fattening. This is shown in Fig. where we compute the static 
potential using c = 0.45, iV = 10 APE-smeared fat links at quenched (3 = 5.85. The loss 
of this part of the potential leads to a suppression of the heavy quark wave function at the 
origin. Although this is an effect that would vanish in the continuum limit (for fixed N, c) , 
it could introduce large scaling violations for short-distance-sensitive quantities. 

To study more directly the effect of fattening on heavy-light decay constants, we have 
computed the decay constants with clover fermions on a 99 lattice subset of quenched set 
CP1, which we call CPF. We have tried four different levels of fattening: c = 0.45 with 



64 



300 — 



> 
0) 



250 



Q 
q— i 



200 



t — r 



"i — i — i — r 



"i — i — i — r 



"i — i — i — r 



"i — i — i — r 



x N f =2 Wilson 
* N f =2 fat clover 



x N f =2 fat clover (corrected) 

const, fit to □ ; CL=0.15 

extrap of □ with quench slope 



¥ .. 













4 



J L 



J I I L 



J I I L 



J I I L 



J I I L 



0.0 0.2 0.4 0.6 

-1 



0.8 



a (GeV) 



FIG. 25: Same as Fig. p^, but for fjj s . 



N = 2, 6, and 10, and c = 0.25 with N = 7. In these cases, csw is set equal to the 
tadpole- improved tree level value 1/wq, with u determined by the plaquette computed with 
the smeared links. The renormalization constants are determined in the same way as for 
the dynamical case (set RF). A comparison of two of the smearing levels with the thin-link 



clover computations is shown in Fig. |2y]. The fat-link fs values are considerably suppressed 
compared to those from the thin links, which in turn are consistent with the results of 
continuum-extrapolated quenched Wilson fermions (see Sec. |TY| ). 

Figure |27| shows that the suppression produced by the lowest and highest levels of fatten- 
ing are consistent. In fact, there is not much difference in the values of the heavy-light decay 
constants among the four different levels of fattening we studied, even though the amount of 
smoothing introduced into the short-distance potential is quite different for the four cases. 
Furthermore, the light-light decay constants with fat clover and thin clover links differ by 
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FIG. 26: Static potential at quenched (3 = 5.85 with and without c = 0.45, N = 10 APE-smeared 
fattening. 



only ~ 7%: Compare the /^-determined lattice spacings of sets CP1 and CPF in Fig. |27|, or 



see Table 0. Note that in the light-light case we are using the correct renormalization fac- 



tors from Ref. JTTJ. This suggests that the ~25% suppression of heavy-light decay constants 
for our fat links may be due more to the use of the incorrect renormalizations (static-light 
instead of heavy-light) than to scaling violations from the smoothing of the short- distance 
potential. Be that as it may, these quenched studies show that the fat clover Nf = 2 results 
may be ignored, at least until fat-link heavy-light renormalization constants are available. 

An alternative approach would be to try to correct the fat-link clover dynamical results 
by the factor (thin link quenched) / (fat link quenched) at a comparable lattice spacing. We 
can do this since the lattice spacings for sets CPF (quenched) and RF (Nf = 2) are quite 
close. (See Table |V].) The corrected fat-link results shown in Figs. [22], |23], and [25] are 
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FIG. 27: Effect of smearing on quenched Jb- The thin clover points are at (3 = 6.0 and 6.15 (sets 
CP1 and J); the fat, at j3 = 6.0 (set CPF). The extrapolation of the thin clover results to the 
continuum is also shown. 

consistent with the Wilson Nf = 2 results. However, we judge that the reasons for the fat- 
link suppression are not well enough understood to be confident that the correction factor 
is the same in the quenched and dynamical cases. We therefore drop the fat clover Nf = 2 
results and use the Wilson results only. 

We emphasize that fat-link actions are formally neither better nor worse than actions with 
thin links — the differences lie in only in the composition and strength of higher dimensional 
(irrelevant) operators. However, from a practical point of view one is interested in actions 
for which particular quantities scale well with lattice spacing. Fat links are intended to 
improve chirality, but chirality is a property of light quarks, not heavy ones. In hindsight, 
there is no physical motivation to construct or use fat link actions for heavy quarks. Some 
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recent developments for fat-link actions for light quarks have been influenced by our 
negative experience — one of the design criteria is to minimize effects such as are shown in 
Fig. ^6|. We are currently studying the behavior of decay constants simulated with thin-link 
heavy quarks and fat-link light quarks. 

C. Partial Quenching and Chiral Extrapolation 

Our central values with Nf = 2 are computed in the "partially quenched" approximation: 
dynamical quark configurations are treated as fixed backgrounds and chiral extrapolation 
is performed in the valence quark mass only. The main justification for using the partially 
quenched approximation can be seen qualitatively in Figs. p2|-p5|: For our range of dynamical 
quark masses and with our statistical and systematic errors, there is no obvious trend in 
the decay constants when the dynamical quark mass is varied at fixed (3. (This statement 
is examined in more detail below.) 

The standard systematic error associated with the valence-mass chiral extrapolation is 
then estimated in exactly the same way as in the quenched approximation (comparison of 
"chiral choice I" with "chiral choice II" — see Sees. |IIID| and[TVl). Effects of chiral logarithms 



at very low quark mass are considered separately in Sec. [VP 



To estimate the systematic error due to partial quenching, we perform a complete addi- 
tional analysis in the "fully unquenched" theory, where the light (u, d) valence quark mass 
on a given lattice set is interpolated or extrapolated to the value of the dynamical mass on 
that set. Since the valence and dynamical quarks are simulated with different lattice actions, 
the equality must be defined by some physical quantity. We demand that the pseudoscalar 
("pion") have the same mass with either action. We then perform chiral extrapolations of 
fq q with m g va i ence = m gjdynamical using data from sets at fixed (3: either (3 = 5.6 (sets G, R, 
S, T, U) or (3 = 5.5 (sets P, M, O, N). Such extrapolations must be performed in physical 
units because they involve different sets with different lattice spacings. To set the scale, we 
use as usual f qq , extrapolated in valence quark mass to the physical u, d point, i.e., f^. Note 
that the scale is set in a partially-quenched manner. However the fully unquenched theory 
is recaptured once the dynamical mass is itself extrapolated to the physical u, d point. 12 



12 This approach could be dangerous if the dependence of f qq on the dynamical quark mass at fixed valence 
mass were so violent that the chiral extrapolation of /q 9 in physical units became uncontrolled. This does 
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In Fig. ^8] we show the chiral extrapolation of fs with m giVa i ence = m 9) d yn amicai at (3 = 5.6. 
We call the dependent variable "/b" because the heavy quark has already been interpolated 
to the b quark mass, as in Sec. |111F| . As independent variable, we use the pseudoscalar 
mass squared, m q , and extrapolate to m q = m 2 . Note that the linear fit is excellent, even 
though it includes very heavy m qq values. However, if we restricted m qq to a safer range for 
a chiral extrapolation (m qq < 0.6(GeV) 2 ), the results would be essentially unchanged. The 
behavior of fo is very similar to that of fs- 
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FIG. 28: "Fully unquenched" chiral extrapolation of fs at j3 = 5.6 (sets G, R, S, T, U). The burst 
shows the extrapolated value when m qq = m 2 . 



Figure £29| shows fs s as a function of the dynamical quark mass at (3 = 5.6. The light 



not appear be the case, as seen in Figs. |2S - 30 below. However, in work in progress |4q| , we employ a safer 
approach, in which the dynamical lattices have matched scales set independently of the valence quarks 
using the static quark potential. 
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valence quark mass has already been interpolated to the strange quark mass, and only the 
dynamical u, d quark mass is varied. With the current statistical and discretization errors, 
there is little evidence here for dynamical quark mass dependence (using an /„. scale.) This 
may be due, at least partially, to staggered flavor violations, which reduce the effective range 
over which the dynamical mass varies. Note, however, that there is a significant difference 
when one compares these dynamical mass points to the infinite mass case (the quenched 



approximation): compare Figs. ^ and |23|. The behavior of fo s is nearly identical to that 
seen in Fig. the other decay constants, such as have similar behavior when they are 
plotted as a function of the dynamical mass for fixed valence mass. 
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FIG. 29: Same as Fig. 28 but for fs s - The valence quark masses do not vary but are held fixed at 
the masses of the b and s. The fit is linear (not a constant), but has quite small slope. 



The chiral extrapolation of /d,//d as a function of dynamical quark mass (represented 
by the dynamical m 2 qq ) is shown in Fig. [30|. For f D , the light valence quark mass is put equal 
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to the dynamical mass; while for fjj B , it is kept equal to the physical strange mass. Since fu a 
has fixed valence quark mass, it, like fs B , changes little with dynamical quark mass; while 
fjj varies more or less linearly, like We therefore fit /d 8 //d to the inverse of a linear 
function in mjL, i.e., to l/(c + dm 2 \ with c and ti allowed to vary The ratio JbJ Jb is fit 
in the same way; while the ratios /b//d s , IbJ Id s , an d /b//d are fit to linear functions. 
(The latter two ratios are, like /b s , almost independent of the dynamical quark mass, and 
so the fitting form makes little difference as long as constant behavior is allowed.) 
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FIG. 30: Chiral extrapolation of /d s //d at /3 = 5.5 (sets P, M, O, N) with the light valence quark 
mass in fjj equal to the dynamical quark mass. The quantity 1/(/d s //d) is fit to a linear function. 



We can now examine the dependence of the fully unquenched quantities on lattice spacing. 
Unfortunately, we can perform the fully unquenched analysis only at two (3 values, 5.5 and 
5.6, for each of which lattice sets exist with four different dynamical quark masses. At the 
third f3 value of our dynamical simulations (/3 = 5.445), we have only a single dynamical 
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mass (am = .025, set L). We attempt a chiral extrapolation there by using the average of 
the (physical) parameters describing the m qq dependence at (3 = 5.6 and 5.5 (as determined 
above). Each parameter has a statistical error estimated by propagating the statistical errors 
of the (3 = 5.6 and 5.5 data, and a systematic error taken to be the difference between the 
average value and the (3 = 5.5 value. The overall error at (3 = 5.445 is then determined 
by adding in quadrature the intrinsic statistical error from set L and the statistical and 
systematic errors coming from the chiral extrapolation. The amount of chiral extrapolation 
required for set L is actually quite small because the physical dynamical quark mass there 
is close to the smallest masses available at (3 — 5.6 and 5.5. Therefore the errors introduced 
by our "synthetic" chiral extrapolation at (3 = 5.445 do not appear to be large. However, 
the fact that the third data point in the fully unquenched analysis must be obtained in this 
way is another reason why we prefer the partially quenched analysis for the central values. 

Figs. [TL| and |32] show the lattice spacing dependence of fs and /b s //b after the fully 
unquenched chiral extrapolations. Like the partially quenched data of Figs. |22| and ^|, the 
fully unquenched data is quite consistent with constant behavior in a. The other decay 
constants and ratios behave similarly. The difference between the result of the constant fits 
in the fully unquenched and partially quenched cases is defined to be the systematic error of 
partial quenching, and is listed for the various quantities in Tables (XI and |X1J|. Given the 



issues in the fully unquenched analysis, we believe that this error determination is merely 
a rough estimate of the magnitude of the effect and do not take the sign of the difference 
seriously. We therefore symmetrize this error in the final error analysis. 

D. Rough Estimate of Chiral Logarithm Effects 

As discussed in Sec. PTEL our rather heavy light mass values preclude a detailed study 



of chiral logarithms. However, an extrapolation of f qq /fQ q (rather than individual decay 
constants), coupled with methods of determining k s and a -1 without significant chiral ex- 
trapolation, should provide an indication of the effect of the logarithms at light quark mass. 
Recall that, in the full theory, the coefficient of the chiral logs in f qq is probably larger than 
in f'Q q . This means that any errors in coming from a quadratic extrapolation of f qq /fQ q 
should be opposite to those in our standard extrapolations of fQ q itself — especially for 
heavy-light decay constant ratios, which are less sensitive to the scale determination. In 
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FIG. 31: Lattice spacing dependence of fs after fully unquenched chiral extrapolation. From left 
to right, the points represent: (3 = 5.6, 5.5, and 5.445. When there is more than one lattice spacing 
at a given /?, the points are plotted at the lattice spacing of the finest lattice (lowest dynamical 
mass). Thus (3 = 5.6 and 5.5 are represented by the lattice spacing of sets R and P, respectively. 
The fit is to a constant. 

particular, the f qq / fQ q approach should overestimate /s s //s, just as our standard approach 
may underestimate it. Indeed, the most significant change from the central value occurs in 
fssl '/b an d is positive. 

Table [VI] of Sec. |I IE] shows the changes in the decay constants with various methods for 
fixing n s and a -1 . Changes in the ratios are given in Table Kill . After eliminating the three 
lines in each table marked with asterisks (see Sec. |I I , we average the changes in decay 
constants and ratios and find the standard deviations of the means. With the exception 
of the quantity fD s /fD, the averages in all cases are positive and larger than the standard 



73 



1.20 



m 



PQ 

«H — I 



1.15 



1.10 



t — i — I — i — i — i — i — I — i — i — i — i — I — i — i — i — i — I — i — i — i — r 



□ N F =2 (fully unquenched) 
const, fit to □; CL=0.98 



J L 



J I I L 



J I I L 



J I I L 



J I I L 



0.0 0.2 0.4 0.6 

-1 



0.8 



a (GeV) 



FIG. 32: Same as Fig. 31, but for f Bs /fB- 



deviations of the means. We define the "error due to chiral logarithm effects" in these cases 
as the signed (positive) number which is the sum of the average and the standard deviation 
of the mean. This is slightly more conservative than just taking the straight average. For 
IdJId, where the average is consistent with 0.00, we take the error as the (unsigned) 
standard deviation of the mean. 

The chiral logarithm effects, while quite significant in the case of Jb s / fs and some of the 
other ratios, appear to be considerably smaller than has been anticipated in Refs. || []]. We 
believe this due to the fact that we set the scale in our central values using /„. and extrapolate 
the light-light and heavy-light decay constants in the same manner. Thus, much of the chiral 
logarithm effects, which are similar in and fs, cancel. 

On the other hand, we emphasize that our estimate of the chiral logarithm effects is, 
for a variety of reasons, rather rough. First of all, the changes in the decay constants and 
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TABLE XIII: Estimates of the effects (in MeV) of chiral logarithms on the extrapolation of decay 
constant ratios. For descriptions of the methods, as well as k s and a -1 , see Table VI in Sec. [IIEj . 
Lines indicated by cl ^ £11*6 eliminated from the averages. 



method 


I'b 3 1 I'b 


Id s 1 Jd 


/b//d s 


fB s /fD 3 


I'b/ Id 


1 


+0.02 


+0.02 


+0.02 


+0.05 


+0.04 




-0.02 


-0.09 


+0.07 


+0.01 


+0.00 


2 


+0.06 


+0.06 


+0.01 


+0.06 


+0.06 




+0.02 


-0.05 


+0.05 


+0.03 


+0.02 


3 


-0.02 


-0.03 


+0.04 


+0.03 


+0.01 




-0.07 


-0.13 


+0.08 


-0.00 


-0.02 


4 


+0.04 


+0.04 


+0.02 


+0.06 


+0.05 




+0.01 


-0.06 


+0.06 


+0.02 


+0.01 


5 


+0.08 


+0.08 


+0.00 


+0.07 


+0.07 




+0.05 


-0.03 


+0.04 


+0.04 


+0.03 


6 


+0.03 


+0.03 


+0.03 


+0.06 


+0.05 




+0.04 


-0.03 


+0.06 


+0.06 


+0.05 


7 


+0.05 


+0.05 


+0.01 


+0.06 


+0.05 




+0.07 


-0.00 


+0.05 


+0.06 


+0.05 


8 


-0.01 


-0.02 


+0.04 


+0.04 


+0.03 




+0.01 


-0.07 


+0.08 


+0.05 


+0.04 


9 


+0.05 


+0.05 


+0.02 


+0.06 


+0.06 




+0.06 


— 0.01 


+0.05 


+0.06 


+0.06 


10 


+0.07 


+0.08 


+0.00 


+0.06 


+0.06 




+0.08 


+0.02 


+0.04 


+0.06 


+0.06 


* 11 


+0.12 


+0.15 


-0.01 


+0.10 


+0.11 


* 


+0.18 


+0.18 


-0.02 


+0.11 


+0.11 


12 


+0.02 


+0.01 


+0.02 


+0.06 


+0.03 


* 


+0.17 


+0.15 


-0.01 


+0.09 


+0.10 


average 


+0.03 


-0.00 


+0.04 


+0.05 


+0.04 


stand, dev. of mean 


0.01 


0.01 


0.01 


0.00 


0.01 
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ratios vary a great deal among the different methods and configurations shown in the tables. 
Indeed, the standard deviation (as opposed to the standard deviation of the mean) of a 
change is typically the same size as the average change and is sometimes larger. Secondly, 
our approach relies on %PT to find the quantities m ss and / ss , and %PT is not necessarily 
rapidly convergent for s quarks. We have also performed only a partially quenched analysis 
of this issue. Because of the size of the errors, we have not attempted to extrapolate the 
dynamical quarks to their physical masses. Finally, we note that there is an inherent (though 
presumably small) inconsistency in our determinations of k s and a -1 , which indirectly use 
the physical values of f„, fx, fn-w and m# (or m v or m^). We cannot force all these quantities 
to have their physical values at once in a theory without a dynamical strange quark. For 
this reason, it is unclear for example whether it is better to use Nf = 2 or Nf = 3 PQ%PT 
in finding / ss ; we hope that our range of methods gives a reasonable range of results. 

Given the crude nature of the chiral log error, we believe that it would be inappropriate 



at this stage to use the computations described in Sec. |I 1 IE| to correct our central values. 
Instead, we use them only for error estimates. 

E. Final Error Estimates and Results 



The magnetic mass error in Tables [X| and |XII] is estimated with almost the same method 



as we used for the quenched calculation. The only difference is that here all the valence 
quarks are of Wilson type, so that there is no reduction of the magnetic mass error in the 
final error budget for the relative number of Wilson and clover estimates. The perturbative 
and 1/M fit errors in the tables are determined in exactly the same manner as in the 
quenched approximation. 

The errors due to finite volume are studied by comparing results on sets R and G, both of 
which have (3 = 5.6 and am = 0.01, but which have spatial volumes 24 3 and 16 3 , respectively. 
Note that all but one of our Nf = 2 sets are large (spatial size ~ 2.1-3.3 fm); only set G 
is comparable in size (~ 1.4 fm) to the quenched lattices. The difference between sets R 
and G is therefore likely to be a considerable overestimate of the actual finite volume error. 
Despite this, the differences are almost never statistically significant. Here a "significant 
difference" is defined as one that is larger than the sum, in quadrature, of the statistical 
errors of the two sets. When the difference is insignificant, we set the finite volume error to 
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zero, as indicated in Tables |X| and |XII| by the notation "— > 0.0." The only case where we 
find a significant (~1.6a) effect is in foj fn- 

The total systematic error within the current approximation (partially quenched Nf = 2 
theory) is then taken to be the sum of all the systematic errors above the line in Tables |XI] and 
pcilj : continuum extrapolation, valence chiral extrapolation, perturbative, magnetic mass, 
1/M fit, and finite volume errors. Since these errors show no evidence of correlations, 
we perform the sum in quadrature. We do, however, treat positive and negative errors 
separately, since the valence chiral extrapolation error represents a binary choice and has a 
well determined sign. 

We still need to estimate the error of the partially quenched Nf = 2 approximation. 
One measure of this error has already been discussed: the partial quenching error. The 
effect of the missing third light virtual quark (the s quark) is estimated in a direct way by 
assuming a simple linear dependence of the decay constants on the number of dynamical 
flavors. The error is thus chosen to be one half the difference of the Nf = 2 and quenched 
calculations. This estimate is labeled "missing dynamical s quark" in Tables |Xj and |XII] . 
We also estimate the effect in two indirect ways: by determining the change in the results 
when (1) the scale is fixed by m p (instead of /„-), and (2) for strange quark quantities, when 
k s is fixed by the vector meson sector (m^) instead of the pseudoscalars. In full QCD (and 
with no other systematic errors!), these differences should vanish, so their size is an estimate 
of the distance we are from the full theory. 

The total error of the partially quenched Nf = 2 approximation is then defined to be 
the maximum of the four estimates below the line in Tables |Xl] and pCIIj : partial quenching, 
scale, k s , and missing dynamical s quark. The latter three estimates have a well-determined 
sign, and we therefore find the maximum positive and maximum negative error separately. 
(As discussed above the partial quenching error is treated symmetrically.) For the individual 
decay constants, the scale and missing dynamical s quark estimates are always largest; while 
the error in the ratios are almost always dominated by the partially quenched error. 

Finally we include an additional error due to the fact that our extrapolations from rather 
large light quark masses cannot see the chiral logarithms directly. This error is estimated 
in Sec. |V H . We emphasize that it is necessarily crude. 

Our final results for heavy-light decay constants, including the effects of dynamical quarks, 
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are: 



f b = 190(7) (±%)(±?)($) MeV; f B , = 217(6)(^)(j*)(U 7 ) MeV 
/ D = 215(6)(tg)(H)(l*) MeV; / Dj = 241(5)(±£) MeV 
^ = 1.16(1)(2)(2)(±*); ^ = 1.14(1)(±1)(3)(1) (46) 

A- = 0.79(2)(U)(3)(^); ^ = 0.92(1) (6) (2) Q 
/A, Jd s 

f -f- = 0.91(2)(±|)(±f)(±g) . 

Here the errors are, respectively, statistical, systematic within the Nf = 2 partially quenched 
approximation, the systematic errors of that approximation (due to partial quenching and 
the missing virtual strange quark), and an estimate of the effect of chiral logarithms. 



The result for fjj a is consistent with experimental results; Ref. |u| obtains f D + = 
280(19)(28)(34) MeV, which the Review of Particle Properties cites as "the best and most 
recent value" f6"4| . Our Nf = 2 values are consistent with recent results of CP-PACS [65 



and preliminary results of JLQCD ||10|| , though our central values of the decay constants 
and ratios /b s //b and fo s l fo are somewhat lower than those of CP-PACS. 

VI. CONCLUSIONS AND FUTURE DIRECTIONS 

Equation |46| and Tables [XI], |X11| , [VT] and [Xlll| summarize our results. Chiral extrapolation, 
continuum extrapolation and perturbation theory are generally the biggest sources of errors 
for the decay constants, while partially quenching, the missing s quark, and the magnetic 
mass are also important for many of the ratios. Because the lattice u, d have necessarily been 
rather heavy, as they have in other lattice calculations to date, the effects of chiral logarithms 
at low quark mass have only been investigated crudely and indirectly. We believe that is 
the error on which we have least control at present. 



Work in progress |45| addresses many of the above issues. Improved actions have de- 
creased the continuum extrapolation errors significantly, as well as eliminated the separate 
magnetic mass error. A dynamical s quark is now explicitly included. Further, since the 
new computations use a wide range of both dynamical and light valence quark masses we 



hope to treat the chiral logarithms explicitly within a partially quenched framework [Bq, 66 



This should provide more direct evidence about the issue about the size of chiral logarithm 
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effects H [|, as well as eliminate the explicit partial quenching error. 

Future calculations will use staggered light quarks, as has already been investigated in 
conjunction with NRQCD heavy quarks |I7||. This will allow for very light valence masses 
and therefore make possible a detailed study of chiral logarithms. To improve the chiral 
extrapolations still more, one-loop chiral perturbation theory calculations that take into 
account staggered taste 13 violation will be needed. Such calculations for pseudoscalar meson 
masses already exist; those for heavy-light decay constants are in progress |68fl . 

The next step after that is likely to involve perturbation theory. Once the other er- 
rors have been reduced, the errors of one-loop perturbative calculations will no longer be 
acceptable. Higher order calculations using automated methods |){J or non-perturbative 
computations will be required. 
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